text:
module ctca
  use, intrinsic :: iso_c_binding
  use mpi
  implicit none

  interface
     function ctca_get_subcomm() bind(c, name='CTCA_get_subcomm')
       import
       integer(kind=4) :: ctca_get_subcomm
     end function ctca_get_subcomm

     function ctca_get_reqinfoitems() bind(c, name='CTCA_get_reqinfoitems')
       import
       integer(kind=4) :: ctca_get_reqinfoitems
     end function ctca_get_reqinfoitems

     subroutine ctcar_init_f() bind(c, name='ctcar_init_f')
       import
     end subroutine ctcar_init_f

     subroutine ctcar_init_detail_f(numareas, numreqs, intparams) bind(c, name='ctcar_init_detail_f')
       import
       integer(kind=4) :: numareas, numreqs, intparams
     end subroutine ctcar_init_detail_f

     subroutine ctcar_regarea_int_f(base, size, areaid) bind(c, name='ctcar_regarea_int_f')
       import
       integer(kind=4) :: areaid
       integer(kind=8) :: size
       integer(kind=4), dimension(1:size) :: base
     end subroutine ctcar_regarea_int_f

     subroutine ctcar_regarea_real4_f(base, size, areaid) bind(c, name='ctcar_regarea_real4_f')
       import
       integer(kind=4) :: areaid
       integer(kind=8) :: size
       real(kind=4), dimension(1:size) :: base
     end subroutine ctcar_regarea_real4_f

     subroutine ctcar_regarea_real8_f(base, size, areaid) bind(c, name='ctcar_regarea_real8_f')
       import
       integer(kind=4) :: areaid
       integer(kind=8) :: size
       real(kind=8), dimension(1:size) :: base
     end subroutine ctcar_regarea_real8_f

     subroutine ctcar_sendreq(intparams, numintparams) bind(c, name='ctcar_sendreq_f')
       import
       integer(kind=4) :: numintparams
       integer(kind=4), dimension(1:numintparams) :: intparams
     end subroutine ctcar_sendreq

     subroutine ctcar_sendreq_hdl(intparams, numintparams, hdl) bind(c, name='ctcar_sendreq_hdl_f')
       import
       integer(kind=4) :: numintparams
       integer(kind=4), dimension(1:numintparams) :: intparams
       integer(kind=8) :: hdl
     end subroutine ctcar_sendreq_hdl

     subroutine ctcar_sendreq_withint_f(intparams, numintparams, data, datanum) &
          bind(c, name='ctcar_sendreq_withint_f')
       import
       integer(kind=4) :: numintparams
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:numintparams) :: intparams
       integer(kind=4), dimension(1:datanum) :: data
     end subroutine ctcar_sendreq_withint_f

     subroutine ctcar_sendreq_withreal4_f(intparams, numintparams, data, datanum) &
          bind(c, name='ctcar_sendreq_withreal4_f')
       import
       integer(kind=4) :: numintparams
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:numintparams) :: intparams
       real(kind=4), dimension(1:datanum) :: data
     end subroutine ctcar_sendreq_withreal4_f

     subroutine ctcar_sendreq_withreal8_f(intparams, numintparams, data, datanum) &
          bind(c, name='ctcar_sendreq_withreal8_f')
       import
       integer(kind=4) :: numintparams
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:numintparams) :: intparams
       real(kind=8), dimension(1:datanum) :: data
     end subroutine ctcar_sendreq_withreal8_f

     subroutine ctcar_sendreq_withint_hdl_f(intparams, numintparams, data, datanum, hdl) &
          bind(c, name='ctcar_sendreq_withint_hdl_f')
       import
       integer(kind=4) :: numintparams
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:numintparams) :: intparams
       integer(kind=4), dimension(1:datanum) :: data
       integer(kind=8) :: hdl
     end subroutine ctcar_sendreq_withint_hdl_f

     subroutine ctcar_sendreq_withreal4_hdl_f(intparams, numintparams, data, datanum, hdl) &
          bind(c, name='ctcar_sendreq_withreal4_hdl_f')
       import
       integer(kind=4) :: numintparams
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:numintparams) :: intparams
       real(kind=4), dimension(1:datanum) :: data
       integer(kind=8) :: hdl
     end subroutine ctcar_sendreq_withreal4_hdl_f

     subroutine ctcar_sendreq_withreal8_hdl_f(intparams, numintparams, data, datanum, hdl) &
          bind(c, name='ctcar_sendreq_withreal8_hdl_f')
       import
       integer(kind=4) :: numintparams
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:numintparams) :: intparams
       real(kind=8), dimension(1:datanum) :: data
       integer(kind=8) :: hdl
     end subroutine ctcar_sendreq_withreal8_hdl_f

     subroutine ctcar_wait(hdl) bind(c, name='ctcar_wait_f')
       import
       integer(kind=8) :: hdl
     end subroutine ctcar_wait

     function ctcar_test_f(hdl) bind(c, name='ctcar_test_f')
       import
       integer(kind=8) :: hdl
       integer(kind=4) :: ctcar_test_f
     end function ctcar_test_f

     subroutine ctcar_finalize() bind(c, name='ctcar_finalize_f')
       import
     end subroutine ctcar_finalize

     subroutine ctcar_prof_start() bind(c, name='ctcar_prof_start_f')
       import
     end subroutine ctcar_prof_start

     subroutine ctcar_prof_stop() bind(c, name='ctcar_prof_stop_f')
       import
     end subroutine ctcar_prof_stop

     subroutine ctcar_prof_start_total() bind(c, name='ctcar_prof_start_total_f')
       import
     end subroutine ctcar_prof_start_total

     subroutine ctcar_prof_stop_total() bind(c, name='ctcar_prof_stop_total_f')
       import
     end subroutine ctcar_prof_stop_total

     subroutine ctcar_prof_start_calc() bind(c, name='ctcar_prof_start_calc_f')
       import
     end subroutine ctcar_prof_start_calc

     subroutine ctcar_prof_stop_calc() bind(c, name='ctcar_prof_stop_calc_f')
       import
     end subroutine ctcar_prof_stop_calc

     subroutine ctcac_init_f() bind(c, name='ctcac_init_f')
       import
     end subroutine ctcac_init_f

     subroutine ctcac_init_detail_f(numareas, numreqs, intparams, bufslotsz, bufslotnum) &
          bind(c, name='ctcac_init_detail_f')
       import
       integer(kind=4) :: numareas, numreqs, intparams, bufslotnum
       integer(kind=8) :: bufslotsz
     end subroutine ctcac_init_detail_f

     subroutine ctcac_regarea_int(areaid) bind(c, name='ctcac_regarea_int_f')
       import
       integer(kind=4) :: areaid
     end subroutine ctcac_regarea_int

     subroutine ctcac_regarea_real4(areaid) bind(c, name='ctcac_regarea_real4_f')
       import
       integer(kind=4) :: areaid
     end subroutine ctcac_regarea_real4

     subroutine ctcac_regarea_real8(areaid) bind(c, name='ctcac_regarea_real8_f')
       import
       integer(kind=4) :: areaid
     end subroutine ctcac_regarea_real8

     function ctcac_isfin_f() bind(c, name='ctcac_isfin_f')
       import
       integer(kind=4) :: ctcac_isfin_f
     end function ctcac_isfin_f

     subroutine ctcac_pollreq(reqinfo, fromrank, intparam, intparamnum) bind(c, name='ctcac_pollreq_f')
       import
       integer(kind=4) :: fromrank, intparamnum
       integer(kind=4), dimension(1:4) :: reqinfo
       integer(kind=4), dimension(1:intparamnum) :: intparam
     end subroutine ctcac_pollreq

     subroutine ctcac_pollreq_withint_f(reqinfo, fromrank, intparam, intparamnum, data, datanum) &
          bind(c, name='ctcac_pollreq_withint_f')
       import
       integer(kind=4) :: fromrank, intparamnum
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:4) :: reqinfo
       integer(kind=4), dimension(1:intparamnum) :: intparam
       integer(kind=4), dimension(1:datanum) :: data
     end subroutine ctcac_pollreq_withint_f

     subroutine ctcac_pollreq_withreal4_f(reqinfo, fromrank, intparam, intparamnum, data, datanum) &
          bind(c, name='ctcac_pollreq_withreal4_f')
       import
       integer(kind=4) :: fromrank, intparamnum
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:4) :: reqinfo
       integer(kind=4), dimension(1:intparamnum) :: intparam
       real(kind=4), dimension(1:datanum) :: data
     end subroutine ctcac_pollreq_withreal4_f

     subroutine ctcac_pollreq_withreal8_f(reqinfo, fromrank, intparam, intparamnum, data, datanum) &
          bind(c, name='ctcac_pollreq_withreal8_f')
       import
       integer(kind=4) :: fromrank, intparamnum
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:4) :: reqinfo
       integer(kind=4), dimension(1:intparamnum) :: intparam
       real(kind=8), dimension(1:datanum) :: data
     end subroutine ctcac_pollreq_withreal8_f

     subroutine ctcac_enqreq(reqinfo, progid, intparam, intparamnum) bind(c, name='ctcac_enqreq_f')
       import
       integer(kind=4) :: progid, intparamnum
       integer(kind=4), dimension(1:4) :: reqinfo
       integer(kind=4), dimension(1:intparamnum) :: intparam
     end subroutine ctcac_enqreq

     subroutine ctcac_enqreq_withint_f(reqinfo, progid, intparam, intparamnum, data, datanum) &
          bind(c, name='ctcac_enqreq_withint_f')
       import
       integer(kind=4) :: progid, intparamnum
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:4) :: reqinfo
       integer(kind=4), dimension(1:intparamnum) :: intparam
       integer(kind=4), dimension(1:datanum) :: data
     end subroutine ctcac_enqreq_withint_f

     subroutine ctcac_enqreq_withreal4_f(reqinfo, progid, intparam, intparamnum, data, datanum) &
          bind(c, name='ctcac_enqreq_withreal4_f')
       import
       integer(kind=4) :: progid, intparamnum
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:4) :: reqinfo
       integer(kind=4), dimension(1:intparamnum) :: intparam
       real(kind=4), dimension(1:datanum) :: data
     end subroutine ctcac_enqreq_withreal4_f

     subroutine ctcac_enqreq_withreal8_f(reqinfo, progid, intparam, intparamnum, data, datanum) &
          bind(c, name='ctcac_enqreq_withreal8_f')
       import
       integer(kind=4) :: progid, intparamnum
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:4) :: reqinfo
       integer(kind=4), dimension(1:intparamnum) :: intparam
       real(kind=8), dimension(1:datanum) :: data
     end subroutine ctcac_enqreq_withreal8_f

     subroutine ctcac_readarea_int_f(areaid, reqrank, offset, size, dest) bind(c, name='ctcac_readarea_int_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       integer(kind=4), dimension(1:size) :: dest
     end subroutine ctcac_readarea_int_f

     subroutine ctcac_readarea_real4_f(areaid, reqrank, offset, size, dest) bind(c, name='ctcac_readarea_real4_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       real(kind=4), dimension(1:size) :: dest
     end subroutine ctcac_readarea_real4_f

     subroutine ctcac_readarea_real8_f(areaid, reqrank, offset, size, dest) bind(c, name='ctcac_readarea_real8_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       real(kind=8), dimension(1:size) :: dest
     end subroutine ctcac_readarea_real8_f

     subroutine ctcac_writearea_int_f(areaid, reqrank, offset, size, src) bind(c, name='ctcac_writearea_int_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       integer(kind=4), dimension(1:size) :: src
     end subroutine ctcac_writearea_int_f

     subroutine ctcac_writearea_real4_f(areaid, reqrank, offset, size, src) bind(c, name='ctcac_writearea_real4_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       real(kind=4), dimension(1:size) :: src
     end subroutine ctcac_writearea_real4_f

     subroutine ctcac_writearea_real8_f(areaid, reqrank, offset, size, src) bind(c, name='ctcac_writearea_real8_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       real(kind=8), dimension(1:size) :: src
     end subroutine ctcac_writearea_real8_f

     subroutine ctcac_finalize() bind(c, name='ctcac_finalize_f')
       import
     end subroutine ctcac_finalize

     subroutine ctcac_prof_start() bind(c, name='ctcac_prof_start_f')
       import
     end subroutine ctcac_prof_start

     subroutine ctcac_prof_stop() bind(c, name='ctcac_prof_stop_f')
       import
     end subroutine ctcac_prof_stop

     subroutine ctcac_prof_start_calc() bind(c, name='ctcac_prof_start_calc_f')
       import
     end subroutine ctcac_prof_start_calc

     subroutine ctcac_prof_stop_calc() bind(c, name='ctcac_prof_stop_calc_f')
       import
     end subroutine ctcac_prof_stop_calc

     subroutine ctcaw_init_f(progid, procspercomm) bind(c, name='ctcaw_init_f')
       import
       integer(kind=4) :: progid, procspercomm
     end subroutine ctcaw_init_f

     subroutine ctcaw_init_detail_f(progid, procspercomm, numareas, intparams) bind(c, name='ctcaw_init_detail_f')
       import
       integer(kind=4) :: progid, procspercomm, numareas, intparams
     end subroutine ctcaw_init_detail_f

     subroutine ctcaw_regarea_int(areaid) bind(c, name='ctcaw_regarea_int_f')
       import
       integer(kind=4) :: areaid
     end subroutine ctcaw_regarea_int

     subroutine ctcaw_regarea_real4(areaid) bind(c, name='ctcaw_regarea_real4_f')
       import
       integer(kind=4) :: areaid
     end subroutine ctcaw_regarea_real4

     subroutine ctcaw_regarea_real8(areaid) bind(c, name='ctcaw_regarea_real8_f')
       import
       integer(kind=4) :: areaid
     end subroutine ctcaw_regarea_real8

     function ctcaw_isfin_f() bind(c, name='ctcaw_isfin_f')
       import
       integer(kind=4) :: ctcaw_isfin_f
     end function ctcaw_isfin_f

     subroutine ctcaw_pollreq(fromrank, intparam, intparamnum) bind(c, name='ctcaw_pollreq_f')
       import
       integer(kind=4) :: fromrank, intparamnum
       integer(kind=4), dimension(1:intparamnum) :: intparam
     end subroutine ctcaw_pollreq

     subroutine ctcaw_pollreq_withint_f(fromrank, intparam, intparamnum, data, datanum) &
          bind(c, name='ctcaw_pollreq_withint_f')
       import
       integer(kind=4) :: fromrank, intparamnum
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:intparamnum) :: intparam
       integer(kind=4), dimension(1:datanum) :: data
     end subroutine ctcaw_pollreq_withint_f

     subroutine ctcaw_pollreq_withreal4_f(fromrank, intparam, intparamnum, data, datanum) &
          bind(c, name='ctcaw_pollreq_withreal4_f')
       import
       integer(kind=4) :: fromrank, intparamnum
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:intparamnum) :: intparam
       real(kind=4), dimension(1:datanum) :: data
     end subroutine ctcaw_pollreq_withreal4_f

     subroutine ctcaw_pollreq_withreal8_f(fromrank, intparam, intparamnum, data, datanum) &
          bind(c, name='ctcaw_pollreq_withreal8_f')
       import
       integer(kind=4) :: fromrank, intparamnum
       integer(kind=8) :: datanum
       integer(kind=4), dimension(1:intparamnum) :: intparam
       real(kind=8), dimension(1:datanum) :: data
     end subroutine ctcaw_pollreq_withreal8_f

     subroutine ctcaw_complete() bind(c, name='ctcaw_complete_f')
       import
     end subroutine ctcaw_complete

     subroutine ctcaw_readarea_int_f(areaid, reqrank, offset, size, dest) bind(c, name='ctcaw_readarea_int_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       integer(kind=4), dimension(1:size) :: dest
     end subroutine ctcaw_readarea_int_f

     subroutine ctcaw_readarea_real4_f(areaid, reqrank, offset, size, dest) bind(c, name='ctcaw_readarea_real4_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       real(kind=4), dimension(1:size) :: dest
     end subroutine ctcaw_readarea_real4_f

     subroutine ctcaw_readarea_real8_f(areaid, reqrank, offset, size, dest) bind(c, name='ctcaw_readarea_real8_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       real(kind=8), dimension(1:size) :: dest
     end subroutine ctcaw_readarea_real8_f

     subroutine ctcaw_writearea_int_f(areaid, reqrank, offset, size, src) bind(c, name='ctcaw_writearea_int_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       integer(kind=4), dimension(1:size) :: src
     end subroutine ctcaw_writearea_int_f

     subroutine ctcaw_writearea_real4_f(areaid, reqrank, offset, size, src) bind(c, name='ctcaw_writearea_real4_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       real(kind=4), dimension(1:size) :: src
     end subroutine ctcaw_writearea_real4_f

     subroutine ctcaw_writearea_real8_f(areaid, reqrank, offset, size, src) bind(c, name='ctcaw_writearea_real8_f')
       import
       integer(kind=4) :: areaid, reqrank
       integer(kind=8) :: offset, size
       real(kind=8), dimension(1:size) :: src
     end subroutine ctcaw_writearea_real8_f

     subroutine ctcaw_finalize() bind(c, name='ctcaw_finalize_f')
       import
     end subroutine ctcaw_finalize

     subroutine ctcaw_prof_start() bind(c, name='ctcaw_prof_start_f')
       import
     end subroutine ctcaw_prof_start

     subroutine ctcaw_prof_stop() bind(c, name='ctcaw_prof_stop_f')
       import
     end subroutine ctcaw_prof_stop

     subroutine ctcaw_prof_start_calc() bind(c, name='ctcaw_prof_start_calc_f')
       import
     end subroutine ctcaw_prof_start_calc

     subroutine ctcaw_prof_stop_calc() bind(c, name='ctcaw_prof_stop_calc_f')
       import
     end subroutine ctcaw_prof_stop_calc

  end interface

  integer(kind=4) :: CTCA_subcomm, CTCAC_REQINFOITEMS
  integer(kind=4), parameter :: CTCA_ROLE_REQ=0
  integer(kind=4), parameter :: CTCA_ROLE_CPL=1
  integer(kind=4), parameter :: CTCA_ROLE_WRK=2

  public CTCA_subcomm, CTCAC_REQINFOITEMS
  public CTCA_ROLE_REQ, CTCA_ROLE_CPL, CTCA_ROLE_WRK
  public CTCAR_init, CTCAR_init_detail
  public CTCAR_regarea_int,  CTCAR_regarea_real4,  CTCAR_regarea_real8
  public CTCAR_sendreq_withint, CTCAR_sendreq_withreal4, CTCAR_sendreq_withreal8
  public CTCAR_sendreq_withint_hdl, CTCAR_sendreq_withreal4_hdl, CTCAR_sendreq_withreal8_hdl
  public CTCAR_test
  public CTCAC_init, CTCAC_init_detail, CTCAC_isfin
  public CTCAC_pollreq_withint, CTCAC_pollreq_withreal4, CTCAC_pollreq_withreal8
  public CTCAC_enqreq_withint, CTCAC_enqreq_withreal4, CTCAC_enqreq_withreal8
  public CTCAC_readarea_int,  CTCAC_readarea_real4,  CTCAC_readarea_real8
  public CTCAW_init, CTCAW_init_detail, CTCAW_isfin
  public CTCAW_pollreq_withint, CTCAW_pollreq_withreal4, CTCAW_pollreq_withreal8
  public CTCAW_readarea_int,  CTCAW_readarea_real4,  CTCAW_readarea_real8

contains
  subroutine ctcar_init()
    implicit none
    call ctcar_init_f()
    CTCA_subcomm = ctca_get_subcomm()
    CTCAC_REQINFOITEMS = ctca_get_reqinfoitems()
  end subroutine ctcar_init

  subroutine ctcar_init_detail(numareas, numreqs, intparams)
    implicit none
    integer(kind=4) :: numareas, numreqs, intparams
    call ctcar_init_detail_f(numareas, numreqs, intparams)
    CTCA_subcomm = ctca_get_subcomm()
    CTCAC_REQINFOITEMS = ctca_get_reqinfoitems()
  end subroutine ctcar_init_detail

  subroutine ctcar_regarea_int(base, size, areaid)
    implicit none
    integer(kind=4) :: areaid
    class(*) :: size
    integer(kind=4), dimension(*) :: base
    integer(kind=8) :: size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
          call ctcar_regarea_int_f(base, size64, areaid)
       type is (integer(kind=8))
          call ctcar_regarea_int_f(base, size, areaid)
       class default
          print *, "Error : ctcar_regqrea_int : wrong data type for size"
    end select
  end subroutine ctcar_regarea_int

  subroutine ctcar_regarea_real4(base, size, areaid)
    implicit none
    integer(kind=4) :: areaid
    class(*) :: size
    real(kind=4), dimension(*) :: base
    integer(kind=8) :: size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
          call ctcar_regarea_real4_f(base, size64, areaid)
       type is (integer(kind=8))
          call ctcar_regarea_real4_f(base, size, areaid)
       class default
          print *, "Error : ctcar_regqrea_real4 : wrong data type for size"
    end select
  end subroutine ctcar_regarea_real4

  subroutine ctcar_regarea_real8(base, size, areaid)
    implicit none
    integer(kind=4) :: areaid
    class(*) :: size
    real(kind=8), dimension(*) :: base
    integer(kind=8) :: size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
          call ctcar_regarea_real8_f(base, size64, areaid)
       type is (integer(kind=8))
          call ctcar_regarea_real8_f(base, size, areaid)
       class default
          print *, "Error : ctcar_regqrea_int : wrong data type for size"
    end select
  end subroutine ctcar_regarea_real8

  subroutine ctcar_sendreq_withint(intparams, numintparams, data, datanum)
    implicit none
    integer(kind=4) :: numintparams
    class(*) :: datanum
    integer(kind=4), dimension(1:numintparams) :: intparams
    integer(kind=4), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcar_sendreq_withint_f(intparams, numintparams, data, datanum64)
       type is (integer(kind=8))
          call ctcar_sendreq_withint_f(intparams, numintparams, data, datanum)
       class default
          print *, "Error : ctcar_sendreq_withint : wrong data type for datanum"
    end select

  end subroutine ctcar_sendreq_withint

  subroutine ctcar_sendreq_withreal4(intparams, numintparams, data, datanum)
    implicit none
    integer(kind=4) :: numintparams
    class(*) :: datanum
    integer(kind=4), dimension(1:numintparams) :: intparams
    real(kind=4), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcar_sendreq_withreal4_f(intparams, numintparams, data, datanum64)
       type is (integer(kind=8))
          call ctcar_sendreq_withreal4_f(intparams, numintparams, data, datanum)
       class default
          print *, "Error : ctcar_sendreq_withreal4 : wrong data type for datanum"
    end select

  end subroutine ctcar_sendreq_withreal4

  subroutine ctcar_sendreq_withreal8(intparams, numintparams, data, datanum)
    implicit none
    integer(kind=4) :: numintparams
    class(*) :: datanum
    integer(kind=4), dimension(1:numintparams) :: intparams
    real(kind=8), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcar_sendreq_withreal8_f(intparams, numintparams, data, datanum64)
       type is (integer(kind=8))
          call ctcar_sendreq_withreal8_f(intparams, numintparams, data, datanum)
       class default
          print *, "Error : ctcar_sendreq_withreal8 : wrong data type for datanum"
    end select

  end subroutine ctcar_sendreq_withreal8

  subroutine ctcar_sendreq_withint_hdl(intparams, numintparams, data, datanum, hdl)
    implicit none
    integer(kind=4) :: numintparams
    class(*) :: datanum
    integer(kind=4), dimension(1:numintparams) :: intparams
    integer(kind=4), dimension(*) :: data
    integer(kind=8) :: datanum64
    integer(kind=8) :: hdl

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcar_sendreq_withint_hdl_f(intparams, numintparams, data, datanum64, hdl)
       type is (integer(kind=8))
          call ctcar_sendreq_withint_hdl_f(intparams, numintparams, data, datanum, hdl)
       class default
          print *, "Error : ctcar_sendreq_withint_hdl : wrong data type for datanum"
    end select

  end subroutine ctcar_sendreq_withint_hdl

  subroutine ctcar_sendreq_withreal4_hdl(intparams, numintparams, data, datanum, hdl)
    implicit none
    integer(kind=4) :: numintparams
    class(*) :: datanum
    integer(kind=4), dimension(1:numintparams) :: intparams
    real(kind=4), dimension(*) :: data
    integer(kind=8) :: datanum64
    integer(kind=8) :: hdl

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcar_sendreq_withreal4_hdl_f(intparams, numintparams, data, datanum64, hdl)
       type is (integer(kind=8))
          call ctcar_sendreq_withreal4_hdl_f(intparams, numintparams, data, datanum, hdl)
       class default
          print *, "Error : ctcar_sendreq_withreal4_hdl : wrong data type for datanum"
    end select

  end subroutine ctcar_sendreq_withreal4_hdl

  subroutine ctcar_sendreq_withreal8_hdl(intparams, numintparams, data, datanum, hdl)
    implicit none
    integer(kind=4) :: numintparams
    class(*) :: datanum
    integer(kind=4), dimension(1:numintparams) :: intparams
    real(kind=8), dimension(*) :: data
    integer(kind=8) :: datanum64
    integer(kind=8) :: hdl

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcar_sendreq_withreal8_hdl_f(intparams, numintparams, data, datanum64, hdl)
       type is (integer(kind=8))
          call ctcar_sendreq_withreal8_hdl_f(intparams, numintparams, data, datanum, hdl)
       class default
          print *, "Error : ctcar_sendreq_withreal8_hdl : wrong data type for datanum"
    end select

  end subroutine ctcar_sendreq_withreal8_hdl

  subroutine ctcac_init()
    implicit none
    call ctcac_init_f()
    CTCA_subcomm = ctca_get_subcomm()
    CTCAC_REQINFOITEMS = ctca_get_reqinfoitems()
  end subroutine ctcac_init

  subroutine ctcac_init_detail(numareas, numreqs, intparams, bufslotsz, bufslotnum)
    implicit none
    integer(kind=4) :: numareas, numreqs, intparams, bufslotnum
    class(*) :: bufslotsz
    integer(kind=8) :: bufslotsz64

    select type (bufslotsz)
       type is (integer(kind=4))
          bufslotsz64 = bufslotsz
          call ctcac_init_detail_f(numareas, numreqs, intparams, bufslotsz64, bufslotnum)
       type is (integer(kind=8))
          call ctcac_init_detail_f(numareas, numreqs, intparams, bufslotsz, bufslotnum)
       class default
          print *, "Error : ctcac_init_detail : wrong data type for bufslotsz"
    end select

    CTCA_subcomm = ctca_get_subcomm()
    CTCAC_REQINFOITEMS = ctca_get_reqinfoitems()
  end subroutine ctcac_init_detail

  function ctcar_test(hdl)
    implicit none
    integer(kind=8) :: hdl
    logical :: ctcar_test

    if (ctcar_test_f(hdl) == 1) then
       ctcar_test = .true.
    else
       ctcar_test = .false.
    end if
  end function ctcar_test

  subroutine ctcac_pollreq_withint(reqinfo, fromrank, intparam, intparamnum, data, datanum)
    implicit none
    integer(kind=4) :: fromrank, intparamnum
    class(*) :: datanum
    integer(kind=4), dimension(1:4) :: reqinfo
    integer(kind=4), dimension(1:intparamnum) :: intparam
    integer(kind=4), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcac_pollreq_withint_f(reqinfo, fromrank, intparam, intparamnum, data, datanum64)
       type is (integer(kind=8))
          call ctcac_pollreq_withint_f(reqinfo, fromrank, intparam, intparamnum, data, datanum)
       class default
          print *, "Error : ctcac_pollreq_withint : wrong data type for datanum"
    end select

  end subroutine ctcac_pollreq_withint

  subroutine ctcac_pollreq_withreal4(reqinfo, fromrank, intparam, intparamnum, data, datanum)
    implicit none
    integer(kind=4) :: fromrank, intparamnum
    class(*) :: datanum
    integer(kind=4), dimension(1:4) :: reqinfo
    integer(kind=4), dimension(1:intparamnum) :: intparam
    real(kind=4), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcac_pollreq_withreal4_f(reqinfo, fromrank, intparam, intparamnum, data, datanum64)
       type is (integer(kind=8))
          call ctcac_pollreq_withreal4_f(reqinfo, fromrank, intparam, intparamnum, data, datanum)
       class default
          print *, "Error : ctcac_pollreq_withreal4 : wrong data type for datanum"
    end select

  end subroutine ctcac_pollreq_withreal4

  subroutine ctcac_pollreq_withreal8(reqinfo, fromrank, intparam, intparamnum, data, datanum)
    implicit none
    integer(kind=4) :: fromrank, intparamnum
    class(*) :: datanum
    integer(kind=4), dimension(1:4) :: reqinfo
    integer(kind=4), dimension(1:intparamnum) :: intparam
    real(kind=8), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcac_pollreq_withreal8_f(reqinfo, fromrank, intparam, intparamnum, data, datanum64)
       type is (integer(kind=8))
          call ctcac_pollreq_withreal8_f(reqinfo, fromrank, intparam, intparamnum, data, datanum)
       class default
          print *, "Error : ctcac_pollreq_withreal8 : wrong data type for datanum"
    end select

  end subroutine ctcac_pollreq_withreal8

  subroutine ctcac_enqreq_withint(reqinfo, progid, intparam, intparamnum, data, datanum)
    implicit none
    integer(kind=4) :: progid, intparamnum
    class(*) :: datanum
    integer(kind=4), dimension(1:4) :: reqinfo
    integer(kind=4), dimension(1:intparamnum) :: intparam
    integer(kind=4), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcac_enqreq_withint_f(reqinfo, progid, intparam, intparamnum, data, datanum64)
       type is (integer(kind=8))
          call ctcac_enqreq_withint_f(reqinfo, progid, intparam, intparamnum, data, datanum)
       class default
          print *, "Error : ctcac_pollreq_withint : wrong data type for datanum"
    end select
  end subroutine ctcac_enqreq_withint

  subroutine ctcac_enqreq_withreal4(reqinfo, progid, intparam, intparamnum, data, datanum)
    implicit none
    integer(kind=4) :: progid, intparamnum
    class(*) :: datanum
    integer(kind=4), dimension(1:4) :: reqinfo
    integer(kind=4), dimension(1:intparamnum) :: intparam
    real(kind=4), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcac_enqreq_withreal4_f(reqinfo, progid, intparam, intparamnum, data, datanum64)
       type is (integer(kind=8))
          call ctcac_enqreq_withreal4_f(reqinfo, progid, intparam, intparamnum, data, datanum)
       class default
          print *, "Error : ctcac_pollreq_withreal4 : wrong data type for datanum"
    end select
  end subroutine ctcac_enqreq_withreal4

  subroutine ctcac_enqreq_withreal8(reqinfo, progid, intparam, intparamnum, data, datanum)
    implicit none
    integer(kind=4) :: progid, intparamnum
    class(*) :: datanum
    integer(kind=4), dimension(1:4) :: reqinfo
    integer(kind=4), dimension(1:intparamnum) :: intparam
    real(kind=8), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcac_enqreq_withreal8_f(reqinfo, progid, intparam, intparamnum, data, datanum64)
       type is (integer(kind=8))
          call ctcac_enqreq_withreal8_f(reqinfo, progid, intparam, intparamnum, data, datanum)
       class default
          print *, "Error : ctcac_pollreq_withreal8 : wrong data type for datanum"
    end select
  end subroutine ctcac_enqreq_withreal8

  subroutine ctcac_readarea_int(areaid, reqrank, offset, size, dest)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    integer(kind=4), dimension(*) :: dest
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcac_readarea_int : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcac_readarea_int : wrong data type for offset"
    end select

    call ctcac_readarea_int_f(areaid, reqrank, offset64, size64, dest)

  end subroutine ctcac_readarea_int

  subroutine ctcac_readarea_real4(areaid, reqrank, offset, size, dest)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    real(kind=4), dimension(*) :: dest
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcac_readarea_real4 : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcac_readarea_real4 : wrong data type for offset"
    end select

    call ctcac_readarea_real4_f(areaid, reqrank, offset64, size64, dest)

  end subroutine ctcac_readarea_real4

  subroutine ctcac_readarea_real8(areaid, reqrank, offset, size, dest)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    real(kind=8), dimension(*) :: dest
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcac_readarea_real8 : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcac_readarea_real8 : wrong data type for offset"
    end select

    call ctcac_readarea_real8_f(areaid, reqrank, offset64, size64, dest)

  end subroutine ctcac_readarea_real8

  subroutine ctcac_writearea_int(areaid, reqrank, offset, size, src)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    integer(kind=4), dimension(*) :: src
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcac_writearea_int : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcac_writearea_int : wrong data type for offset"
    end select

    call ctcac_writearea_int_f(areaid, reqrank, offset64, size64, src)

  end subroutine ctcac_writearea_int

  subroutine ctcac_writearea_real4(areaid, reqrank, offset, size, src)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    real(kind=4), dimension(*) :: src
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcac_writearea_real4 : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcac_writearea_real4 : wrong data type for offset"
    end select

    call ctcac_writearea_real4_f(areaid, reqrank, offset64, size64, src)

  end subroutine ctcac_writearea_real4

  subroutine ctcac_writearea_real8(areaid, reqrank, offset, size, src)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    real(kind=8), dimension(*) :: src
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcac_writearea_real8 : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcac_writearea_real8 : wrong data type for offset"
    end select

    call ctcac_writearea_real8_f(areaid, reqrank, offset64, size64, src)

  end subroutine ctcac_writearea_real8

  function ctcac_isfin()
    implicit none
    logical :: ctcac_isfin

    if (ctcac_isfin_f() == 1) then
       ctcac_isfin = .true.
    else
       ctcac_isfin = .false.
    end if
  end function ctcac_isfin

  subroutine ctcaw_init(progid, procspercomm)
    implicit none
    integer(kind=4) :: progid, procspercomm
    call ctcaw_init_f(progid, procspercomm)
    CTCA_subcomm = ctca_get_subcomm()
    CTCAC_REQINFOITEMS = ctca_get_reqinfoitems()
  end subroutine ctcaw_init

  subroutine ctcaw_init_detail(progid, procspercomm, numareas, intparams)
    implicit none
    integer(kind=4) :: progid, procspercomm, numareas, intparams
    call ctcaw_init_detail_f(progid, procspercomm, numareas, intparams)
    CTCA_subcomm = ctca_get_subcomm()
    CTCAC_REQINFOITEMS = ctca_get_reqinfoitems()
  end subroutine ctcaw_init_detail

  subroutine ctcaw_pollreq_withint(fromrank, intparam, intparamnum, data, datanum)
    implicit none
    integer(kind=4) :: fromrank, intparamnum
    class(*) :: datanum
    integer(kind=4), dimension(1:intparamnum) :: intparam
    integer(kind=4), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcaw_pollreq_withint_f(fromrank, intparam, intparamnum, data, datanum64)
       type is (integer(kind=8))
          call ctcaw_pollreq_withint_f(fromrank, intparam, intparamnum, data, datanum)
       class default
          print *, "Error : ctcac_pollreq_withint : wrong data type for datanum"
    end select
  end subroutine ctcaw_pollreq_withint

  subroutine ctcaw_pollreq_withreal4(fromrank, intparam, intparamnum, data, datanum)
    implicit none
    integer(kind=4) :: fromrank, intparamnum
    class(*) :: datanum
    integer(kind=4), dimension(1:intparamnum) :: intparam
    real(kind=4), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcaw_pollreq_withreal4_f(fromrank, intparam, intparamnum, data, datanum64)
       type is (integer(kind=8))
          call ctcaw_pollreq_withreal4_f(fromrank, intparam, intparamnum, data, datanum)
       class default
          print *, "Error : ctcac_pollreq_withreal4 : wrong data type for datanum"
    end select
  end subroutine ctcaw_pollreq_withreal4

  subroutine ctcaw_pollreq_withreal8(fromrank, intparam, intparamnum, data, datanum)
    implicit none
    integer(kind=4) :: fromrank, intparamnum
    class(*) :: datanum
    integer(kind=4), dimension(1:intparamnum) :: intparam
    real(kind=8), dimension(*) :: data
    integer(kind=8) :: datanum64

    select type (datanum)
       type is (integer(kind=4))
          datanum64 = datanum
          call ctcaw_pollreq_withreal8_f(fromrank, intparam, intparamnum, data, datanum64)
       type is (integer(kind=8))
          call ctcaw_pollreq_withreal8_f(fromrank, intparam, intparamnum, data, datanum)
       class default
          print *, "Error : ctcac_pollreq_withreal8 : wrong data type for datanum"
    end select
  end subroutine ctcaw_pollreq_withreal8

  subroutine ctcaw_readarea_int(areaid, reqrank, offset, size, dest)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    integer(kind=4), dimension(*) :: dest
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcaw_readarea_int : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcaw_readarea_int : wrong data type for offset"
    end select

    call ctcaw_readarea_int_f(areaid, reqrank, offset64, size64, dest)

  end subroutine ctcaw_readarea_int

  subroutine ctcaw_readarea_real4(areaid, reqrank, offset, size, dest)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    real(kind=4), dimension(*) :: dest
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcaw_readarea_real4 : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcaw_readarea_real4 : wrong data type for offset"
    end select

    call ctcaw_readarea_real4_f(areaid, reqrank, offset64, size64, dest)

  end subroutine ctcaw_readarea_real4

  subroutine ctcaw_readarea_real8(areaid, reqrank, offset, size, dest)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    real(kind=8), dimension(*) :: dest
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcaw_readarea_real8 : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcaw_readarea_real8 : wrong data type for offset"
    end select

    call ctcaw_readarea_real8_f(areaid, reqrank, offset64, size64, dest)

  end subroutine ctcaw_readarea_real8

  subroutine ctcaw_writearea_int(areaid, reqrank, offset, size, src)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    integer(kind=4), dimension(*) :: src
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcaw_writearea_int : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcaw_writearea_int : wrong data type for offset"
    end select

    call ctcaw_writearea_int_f(areaid, reqrank, offset64, size64, src)

  end subroutine ctcaw_writearea_int

  subroutine ctcaw_writearea_real4(areaid, reqrank, offset, size, src)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    real(kind=4), dimension(*) :: src
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcaw_writearea_real4 : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcaw_writearea_real4 : wrong data type for offset"
    end select

    call ctcaw_writearea_real4_f(areaid, reqrank, offset64, size64, src)

  end subroutine ctcaw_writearea_real4

  subroutine ctcaw_writearea_real8(areaid, reqrank, offset, size, src)
    implicit none
    integer(kind=4) :: areaid, reqrank
    class(*) :: offset, size
    real(kind=8), dimension(*) :: src
    integer(kind=8) :: offset64, size64

    select type (size)
       type is (integer(kind=4))
          size64 = size
       type is (integer(kind=8))
          size64 = size
       class default
          print *, "Error : ctcaw_writearea_real8 : wrong data type for size"
    end select

    select type (offset)
       type is (integer(kind=4))
          offset64 = offset
       type is (integer(kind=8))
          offset64 = offset
       class default
          print *, "Error : ctcaw_writearea_real8 : wrong data type for offset"
    end select

    call ctcaw_writearea_real8_f(areaid, reqrank, offset64, size64, src)

  end subroutine ctcaw_writearea_real8

  function ctcaw_isfin()
    implicit none
    logical :: ctcaw_isfin

    if (ctcaw_isfin_f() == 1) then
       ctcaw_isfin = .true.
    else
       ctcaw_isfin = .false.
    end if
  end function ctcaw_isfin

end module ctca

max_stars_repo_name: tnanri/cotocoa

max_issues_repo_name: tnanri/cotocoa

max_forks_repo_name: tnanri/cotocoa


text:
      program t
      goto (100,200) 1
      print *,'bad'
 100  print *,100
      goto 999
 200  print *,200
 999  print *,'done'
      end program t

ext: f

max_stars_repo_name: maddenp/ppp

max_issues_repo_name: maddenp/ppp

max_forks_repo_name: maddenp/ppp


! PR middle-end/82004
! { dg-do run }
! { dg-options "-Ofast" }

  integer, parameter :: r8 = selected_real_kind(13), i4 = kind(1)
  integer (i4), parameter :: a = 400, b = 2
  real (r8), parameter, dimension(b) :: c = (/ .001_r8, 10.00_r8 /)
  real (r8) :: d, e, f, g, h
  real (r8), parameter :: j &
    = 10**(log10(c(1))-(log10(c(b))-log10(c(1)))/real(a))

  d = c(1)
  e = c(b)
  f = (log10(e)-log10(d))/real(a)
  g = log10(d) - f
  h = 10**(g)
  if (h.ne.j) stop 1
end


! { dg-require-effective-target tls }
      module omp_threadprivate1
	common /T/ a
      end module omp_threadprivate1
      subroutine bad1
	use omp_threadprivate1
!$omp threadprivate (/T/)	! { dg-error "not found" }
      end subroutine bad1
      subroutine bad2
	common /S/ b
!$omp threadprivate (/S/)
      contains
	subroutine bad3
!$omp parallel copyin (/T/)	! { dg-error "not found" }
!$omp end parallel		! { dg-error "" }
	end subroutine bad3
      end subroutine bad2


text:
      program Haghighat

      integer a, b, apb, amb

      read *, a, b
      apb = a + b
      amb = a - b
      read *, b
      print *, apb, (apb + amb)/2
      end


text:
module pksdft_fft
use constants, only: pi
use types, only: dp
use utils, only: assert, stop_error, clock, allocate_mold
use pofdft_fft, only: preal2fourier, pfourier2real
use arpack, only: peig
implicit none
private
public solve_schroedinger

contains

subroutine solve_schroedinger(myid, comm_all, commy, commz, Ng, nsub, Vloc, &
        L, G2, cutfn, nev, ncv, eigs, orbitals)
integer, intent(in) :: myid, comm_all, commy, commz, Ng(3), nsub(3), nev, ncv
real(dp), intent(in) :: Vloc(:,:,:) ! Local effective potential
real(dp), intent(in) :: G2(:,:,:), L(:), cutfn(:,:,:)
real(dp), intent(out) :: eigs(:) ! eigs(nev)
! orbitals(Ng_local(1),Ng_local(2),Ng_local(3),nev)
real(dp), intent(out) :: orbitals(:,:,:,:)

real(dp), allocatable :: d(:), v(:,:)
integer :: Ng_local(3), n, i
real(dp) :: t1, t2
logical :: verbose
verbose = .false.
Ng_local = Ng / nsub
n = product(Ng_local)
allocate(v(n,ncv), d(ncv))
call cpu_time(t1)
call peig(comm_all, myid, n, nev, ncv, "SA", av, d, v)
call cpu_time(t2)
eigs = d(:nev)+10
orbitals = reshape(v(:,:nev), [Ng_local(1),Ng_local(2),Ng_local(3),nev]) &
    * sqrt(product(Ng/L))
if (myid == 0 .and. verbose) then
    print *, "Arpack Time:", t2-t1
end if

contains

    subroutine av(x, y)
    ! Compute y = A*x
    real(dp), intent(in) :: x(:)
    real(dp), intent(out) :: y(:)
    complex(dp), dimension(Ng_local(1),Ng_local(2),Ng_local(3)) :: &
        psi, psiG, psiG_vloc
    call preal2fourier(reshape(x, [Ng_local(1),Ng_local(2),Ng_local(3)]), &
        psiG, commy, commz, Ng, nsub)
    psiG = psiG * cutfn
    ! psiG is our starting point

    ! Apply kinetic and potential
    call pfourier2real(psiG, psi, commy, commz, Ng, nsub)
    call preal2fourier((Vloc-10)*psi, psiG_vloc, commy, commz, Ng, nsub)
    psiG_vloc = psiG_vloc * cutfn

    psiG = G2/2*psiG + psiG_vloc

    ! Convert to real space at the end
    call pfourier2real(psiG, psi, commy, commz, Ng, nsub)
    y = reshape(real(psi,dp), [product(Ng_local)])
    end

end subroutine

end module

max_stars_repo_name: certik/hfsolver

max_issues_repo_name: certik/hfsolver

max_forks_repo_name: certik/hfsolver


program HDF5_shapes
!! This program shows how HDF5 dimension orders are distinct in different langauges
use hdf5_interface, only: hdf5_file,hsize_t
use, intrinsic:: iso_fortran_env, only: real64, stdout=>output_unit
implicit none

type(hdf5_file) :: h5f
character(1024) :: argv
character(:), allocatable :: fn, dname, cmd
integer(HSIZE_T), allocatable :: dims(:)
integer :: ierr

if (command_argument_count() /= 2) error stop "filename dset_name"

call get_command_argument(1, argv)
fn = trim(argv)

call get_command_argument(2, argv)
dname = trim(argv)

call h5f%initialize(fn, status='old', action='r')

call h5f%shape(dname, dims)

print '(/,A,100I8)', 'Fortran dims: ',dims

call h5F%finalize()

!> use HDF5 tool to see their view
write(stdout,'(/,A)', advance='no') 'H5LS: '
call execute_command_line('h5ls '//fn//' | grep '//dname)

!> Python
call execute_command_line('which python > /dev/null', exitstat=ierr)
if (ierr==0) then
  write(stdout,'(/,A)', advance='no') 'Python h5py: '
   cmd = "import h5py; f = h5py.File('"//fn//"','r'); print(f['"//dname//"'].shape)"
  call execute_command_line('python -c "'//cmd//'"', exitstat=ierr)
endif
if (ierr/=0) print *,'Python HDF5 not available. '

!> Matlab
call execute_command_line('which matlab > /dev/null', exitstat=ierr)
if (ierr==0) then
  write(stdout,'(/,A)', advance='no') 'Matlab dims: '
   cmd="i=h5info('"//fn//"','/"//dname//"'); disp(i.Dataspace.Size); pause(0.1); exit"
  call execute_command_line('matlab -nojvm -r "'//cmd//'" | tail -n3', wait=.true., exitstat=ierr)
endif
if (ierr/=0) print *,'Matlab HDF5 not available. '

end program


text:
*     TB01ID EXAMPLE PROGRAM TEXT.
*     Copyright (c) 2002-2020 NICONET e.V.
*
*     .. Parameters ..
      INTEGER          NIN, NOUT
      PARAMETER        ( NIN = 5, NOUT = 6 )
      INTEGER          NMAX, MMAX, PMAX
      PARAMETER        ( NMAX = 20, MMAX = 20, PMAX = 20 )
      INTEGER          LDA, LDB, LDC
      PARAMETER        ( LDA = NMAX, LDB = NMAX, LDC = PMAX )
*     .. Local Scalars ..
      CHARACTER*1      JOB
      INTEGER          I, INFO, J, M, N, P
      DOUBLE PRECISION MAXRED
*     .. Local Arrays ..
      DOUBLE PRECISION A(LDA,NMAX), B(LDB,MMAX), C(LDC,NMAX),
     $                 SCALE(NMAX)
*     .. External Subroutines ..
      EXTERNAL         TB01ID, UD01MD
*     .. Executable Statements ..
*
      WRITE ( NOUT, FMT = 99999 )
*     Skip the heading in the data file and read the data.
      READ ( NIN, FMT = '()' )
      READ ( NIN, FMT = * ) N, M, P, JOB, MAXRED
      IF ( N.LT.0 .OR. N.GT.NMAX ) THEN
         WRITE ( NOUT, FMT = 99993 ) N
      ELSE
         READ ( NIN, FMT = * ) ( ( A(I,J), J = 1,N ), I = 1,N )
         IF ( M.LT.0 .OR. M.GT.MMAX ) THEN
            WRITE ( NOUT, FMT = 99992 ) M
         ELSE
            READ ( NIN, FMT = * ) ( ( B(I,J), J = 1,M ), I = 1,N )
            IF ( P.LT.0 .OR. P.GT.MMAX ) THEN
               WRITE ( NOUT, FMT = 99991 ) P
            ELSE
               READ ( NIN, FMT = * ) ( ( C(I,J), J = 1,N ), I = 1,P )
*              Balance system matrix S.
               CALL TB01ID( JOB, N, M, P, MAXRED, A, LDA, B, LDB, C,
     $                      LDC, SCALE, INFO )
*
               IF ( INFO.NE.0 ) THEN
                  WRITE ( NOUT, FMT = 99998 ) INFO
               ELSE
                  CALL UD01MD( N, N, 5, NOUT, A, LDA,
     $                        'The balanced matrix A', INFO )
                  IF ( M.GT.0 )
     $               CALL UD01MD( N, M, 5, NOUT, B, LDB,
     $                            'The balanced matrix B', INFO )
                  IF ( P.GT.0 )
     $               CALL UD01MD( P, N, 5, NOUT, C, LDC,
     $                            'The balanced matrix C', INFO )
                  CALL UD01MD( 1, N, 5, NOUT, SCALE, 1,
     $                        'The scaling vector SCALE', INFO )
                  WRITE ( NOUT, FMT = 99994 ) MAXRED
               END IF
            END IF
         END IF
      END IF
      STOP
*
99999 FORMAT (' TB01ID EXAMPLE PROGRAM RESULTS',/1X)
99998 FORMAT (' INFO on exit from TB01ID = ',I2)
99994 FORMAT (/' MAXRED is ',E13.4)
99993 FORMAT (/' N is out of range.',/' N = ',I5)
99992 FORMAT (/' M is out of range.',/' M = ',I5)
99991 FORMAT (/' P is out of range.',/' P = ',I5)
      END

max_stars_repo_name: bnavigator/SLICOT-Reference

max_issues_repo_name: bnavigator/SLICOT-Reference

max_forks_repo_name: bnavigator/SLICOT-Reference


text:
!----------------------------------------------------------------
! (c) Copyright, 2016 by the Regents of the University of California.
! EMfldclass: ElectroMagnetic field data container class
!             in Lattice module of APPLICATION layer.
! Version: 1.0
! Author: Ji Qiang, LBNL
! Description: This class contains discrete EM field data (as a function of
!              x,y,z) or (r,z) and analytical representation of EM field data (user can
!              supply the function form). The linear transfer map is also
!              computed base on the field on the axis.
! Comments:
!----------------------------------------------------------------
      module EMfldclass
        use PhysConstclass
        use Dataclass
        integer, private, parameter :: Nparam = 14
        type EMfld
          !Itype = 110
          integer :: Nseg,Mapstp,Itype
          double precision :: Length
          double precision, dimension(Nparam) :: Param
          ! Param(1) : zedge
          !      (2) : scale
          !      (3) : RF frequency
          !      (4) : theta0
          !      (5) : file ID
          !      (6) : x radius
          !      (7) : y radius
          !      (8) : x misalignment error
          !      (9) : y misalignment error
          !      (10) : rotation error x
          !      (11) : rotation error y
          !      (12) : rotation error z
          !      (13) : flag for 3D discrete data(1),analytical+discrete(2),analytical(other)
          !      (14) : flag for Cartisian(2) or Cylindrical coordintate(1)
        end type EMfld
        interface getparam_EMfld
          module procedure getparam1_EMfld,  &
                          getparam2_EMfld,   &
                          getparam3_EMfld
        end interface
        interface setparam_EMfld
          module procedure setparam1_EMfld,  &
                           setparam2_EMfld, setparam3_EMfld
        end interface
      contains
        subroutine construct_EMfld(this,numseg,nmpstp,type,blength)
        implicit none
        type (EMfld), intent(out) :: this
        integer, intent(in) :: numseg,nmpstp,type
        double precision, intent(in) :: blength

        this%Nseg = numseg
        this%Mapstp = nmpstp
        this%Itype = type
        this%Length = blength

        this%Param = 0.0

        end subroutine construct_EMfld

        subroutine setparam1_EMfld(this,i,value)
        implicit none
        type (EMfld), intent(inout) :: this
        integer, intent(in) :: i
        double precision, intent(in) :: value

        this%Param(i) = value

        end subroutine setparam1_EMfld

        subroutine setparam2_EMfld(this,values)
        implicit none
        type (EMfld), intent(inout) :: this
        double precision, dimension(:), intent(in) :: values

        this%Param = values

        end subroutine setparam2_EMfld

        subroutine setparam3_EMfld(this,numseg,nmpstp,type,blength)
        implicit none
        type (EMfld), intent(inout) :: this
        integer, intent(in) :: numseg,nmpstp,type
        double precision, intent(in) :: blength

        this%Nseg = numseg
        this%Mapstp = nmpstp
        this%Itype = type
        this%Length = blength

        end subroutine setparam3_EMfld

        subroutine getparam1_EMfld(this,i,blparam)
        implicit none
        type (EMfld), intent(in) :: this
        integer, intent(in) :: i
        double precision, intent(out) :: blparam

        blparam = this%Param(i)

        end subroutine getparam1_EMfld

        subroutine getparam2_EMfld(this,blparams)
        implicit none
        type (EMfld), intent(in) :: this
        double precision, dimension(:), intent(out) :: blparams

        blparams = this%Param

        end subroutine getparam2_EMfld

        subroutine getparam3_EMfld(this,blength,bnseg,bmapstp,&
                                       btype)
        implicit none
        type (EMfld), intent(in) :: this
        double precision, intent(out) :: blength
        integer, intent(out) :: bnseg,bmapstp,btype

        blength = this%Length
        bnseg = this%Nseg
        bmapstp = this%Mapstp
        btype = this%Itype

        end subroutine getparam3_EMfld

        subroutine maplinear_EMfld(t,tau,xm,this,refpt,Bchg,Bmass)
        implicit none
        include 'mpif.h'
        double precision, intent(in) :: t,tau,Bchg,Bmass
        double precision, dimension(6,6), intent(out) :: xm
        double precision, dimension(6), intent(inout) :: refpt
        type (EMfld), intent(in) :: this
        double precision, dimension(18) :: y
        double precision :: h,t0,gi,betai,ui,squi,sqi3
        double precision :: dlti,thli,gf,betaf,uf,squf,sqf3
        double precision :: dltf,thlf,zedge,escale,ww,theta0,qmcc
        double precision :: ein,e1in,e2in,ef,e1f,e2f
        double precision :: sinthi,costhi,uprimi,qpwi,tfin,sinthf,costhf
        double precision :: uprimf,qpwf
        integer :: mpstp

        xm = 0.0

        zedge = this%Param(1)
        escale = this%Param(2)
        ww = this%Param(3)/Scfreq
        theta0 = this%Param(4)*asin(1.0)/90
        qmcc = Bchg/Bmass
        mpstp = this%Mapstp

        y(1)=0.0
        y(2)=0.0
        y(3)=0.0
        y(4)=0.0
        y(5)=refpt(5)
        y(6)=refpt(6)
        y(7)=1.0
        y(8)=0.0
        y(9)=0.0
        y(10)=1.0
        y(11)=1.0
        y(12)=0.0
        y(13)=0.0
        y(14)=1.0
        y(15)=1.0
        y(16)=0.0
        y(17)=0.0
        y(18)=1.0
        h=tau/mpstp
        t0=t
        gi=-refpt(6)
        betai=sqrt(gi**2-1.)/gi
        ui=gi*betai
        squi=sqrt(ui)
        sqi3=sqrt(ui)**3

!        call getaxfldE_EMfld(t,this,ein,e1in,e2in)
        call getaxfldEfc_EMfld(t,this,ein,e1in,e2in)
        sinthi=sin(ww*refpt(5)+theta0)
        costhi=cos(ww*refpt(5)+theta0)
        uprimi=qmcc*ein/betai*costhi
        qpwi=0.5*qmcc*Scxl/(ui*ww)
        dlti=Scxl*(0.5*uprimi/ui-qpwi*e1in*sinthi)
        thli=1.5*Scxl*(uprimi/ui)

        call rk6i_EMfld(h,mpstp,t0,y,18,this,Bchg,Bmass)

        refpt(5)=y(5)
        refpt(6)=y(6)
        gf=-refpt(6)
        betaf=sqrt(gf**2-1.)/gf
        uf=gf*betaf
        squf=sqrt(uf)
        sqf3=sqrt(uf)**3
        tfin = t + tau
!        call getaxfldE_EMfld(tfin,this,ef,e1f,e2f)
        call getaxfldEfc_EMfld(tfin,this,ef,e1f,e2f)
        sinthf=sin(ww*refpt(5)+theta0)
        costhf=cos(ww*refpt(5)+theta0)
        uprimf=qmcc*ef/betaf*costhf
        qpwf=0.5*qmcc*Scxl/(uf*ww)
        dltf=Scxl*(0.5*uprimf/uf-qpwf*e1f*sinthf)
        thlf=1.5*Scxl*(uprimf/uf)

        xm(1,1)= y(7)*squi/squf+y(9)*squi/squf*dlti
        xm(2,1)=(y(8)-y(7)*dltf)*squi*squf+(y(10)-y(9)*dltf)*squi*squf*&
                 dlti
        xm(1,2)= y(9)/(squi*squf)
        xm(2,2)=(y(10)-y(9)*dltf)*squf/squi
        xm(3,3)= y(11)*squi/squf+y(13)*squi/squf*dlti
        xm(4,3)= &
        (y(12)-y(11)*dltf)*squi*squf+(y(14)-y(13)*dltf)*squi*squf*dlti
        xm(3,4)= y(13)/(squi*squf)
        xm(4,4)=(y(14)-y(13)*dltf)*squf/squi
        xm(5,5)= y(15)*sqi3/sqf3+y(17)*sqi3/sqf3*thli
        xm(6,5)= &
        (y(16)-y(15)*thlf)*sqi3*sqf3+(y(18)-y(17)*thlf)*sqi3*sqf3*thli
        xm(5,6)= y(17)/(sqi3*sqf3)
        xm(6,6)=(y(18)-y(17)*thlf)*sqf3/sqi3

        end subroutine maplinear_EMfld

        subroutine rk6i_EMfld(h,ns,t,y,nvar,this,Bchg,Bmass)
        implicit none
        include 'mpif.h'
        integer, intent(in) :: ns,nvar
        double precision, intent(inout) :: t
        double precision, intent(in) :: h,Bchg,Bmass
        double precision, dimension(nvar), intent(inout) :: y
        type (EMfld), intent(in) :: this
        double precision, dimension(nvar) :: yt,a,b,c,d,e,f,g,o,p
        double precision:: tint,tt
        integer :: i
        tint=t
        do i=1,ns
          call intfunc1_EMfld(t,y,f,this,Bchg,Bmass)
          a=h*f
          yt=y+a/9.d0
          tt=t+h/9.d0
          call intfunc1_EMfld(tt,yt,f,this,Bchg,Bmass)
          b=h*f
          yt=y + (a + 3.d0*b)/24.d0
          tt=t+h/6.d0
          call intfunc1_EMfld(tt,yt,f,this,Bchg,Bmass)
          c=h*f
          yt=y+(a-3.d0*b+4.d0*c)/6.d0
          tt=t+h/3.d0
          call intfunc1_EMfld(tt,yt,f,this,Bchg,Bmass)
          d=h*f
          yt=y + (278.d0*a - 945.d0*b + 840.d0*c + 99.d0*d)/544.d0
          tt=t+.5d0*h
          call intfunc1_EMfld(tt,yt,f,this,Bchg,Bmass)
          e=h*f
          yt=y + (-106.d0*a+273.d0*b-104.d0*c-107.d0*d+48.d0*e)/6.d0
          tt = t+2.d0*h/3.d0
          call intfunc1_EMfld(tt,yt,f,this,Bchg,Bmass)
          g=h*f
          yt = y+(110974.d0*a-236799.d0*b+68376.d0*c+ &
                  103803.d0*d-10240.d0*e + 1926.d0*g)/45648.d0
          tt = t + 5.d0*h/6.d0
          call intfunc1_EMfld(tt,yt,f,this,Bchg,Bmass)
          o=h*f
          yt = y+(-101195.d0*a+222534.d0*b-71988.d0*c-  &
                  26109.d0*d-2.d4*e-72.d0*g+22824.d0*o)/25994.d0
          tt = t + h
          call intfunc1_EMfld(tt,yt,f,this,Bchg,Bmass)
          p=h*f
          y = y+(41.d0*a+216.d0*c+27.d0*d+ &
                 272.d0*e+27.d0*g+216.d0*o+41.d0*p)/840.d0
          t=tint+i*h
        enddo

        end subroutine rk6i_EMfld

        subroutine intfunc1_EMfld(t,y,f,this,Bchg,Bmass)
        implicit none
        include 'mpif.h'
        double precision, intent(in) :: t,Bchg,Bmass
        double precision, dimension(:), intent(in) :: y
        double precision, dimension(:), intent(out) :: f
        type (EMfld), intent(in) :: this
        double precision :: gamma0,beta0,gbet,s11,s33,s55,s11tmp
        double precision :: zedge,escale,ez1,ezp1,ezpp1,ww,theta0,qmcc,&
                            sinphi,cosphi,rfdsgn
        integer :: my_rank, ierr

        zedge = this%Param(1)
!        call getaxfldE_EMfld(t,this,ez1,ezp1,ezpp1)
        call getaxfldEfc_EMfld(t,this,ez1,ezp1,ezpp1)
        ww = this%Param(3)/Scfreq
        theta0 = this%Param(4)*asin(1.0)/90
        qmcc = Bchg/Bmass

        f(1) = 0.0
        f(2) = 0.0
        f(3) = 0.0
        f(4) = 0.0

        ! synchronous particle:
        gamma0=-y(6)
        beta0=sqrt(gamma0**2-1.)/gamma0
        gbet=sqrt((gamma0-1.0)*(gamma0+1.0))
        sinphi=sin(ww*y(5)+theta0)
        cosphi=cos(ww*y(5)+theta0)
        rfdsgn=ez1*cosphi
        f(5)=1.0/(beta0*Scxl)
        f(6)=-qmcc*rfdsgn

        ! matrix elements
        s11tmp=0.5e0*(1.0+0.5*gamma0**2)* &
         (qmcc*rfdsgn/beta0**2/gamma0**2)**2 + &
          qmcc*0.5/beta0**3/gamma0**3*ez1*sinphi*ww/Scxl
        s11=s11tmp*Scxl
        s33=s11tmp*Scxl
        s55=  &
         -1.5e0*qmcc/beta0**2/gamma0*ezp1*cosphi &
         +(beta0**2+0.5)/beta0**3/gamma0*qmcc*ez1*sinphi*ww/Scxl &
         +1.5e0*(1.0-0.5*gamma0**2)* &
         (qmcc*rfdsgn/beta0**2/gamma0**2)**2
        s55=Scxl*s55

        f(7)=y(8)/Scxl
        f(8)=-s11*y(7)
        f(9)=y(10)/Scxl
        f(10)=-s11*y(9)
        f(11)=y(12)/Scxl
        f(12)=-s33*y(11)
        f(13)=y(14)/Scxl
        f(14)=-s33*y(13)
        f(15)=y(16)/Scxl
        f(16)=-s55*y(15)
        f(17)=y(18)/Scxl
        f(18)=-s55*y(17)

        end subroutine intfunc1_EMfld

        !interpolate the field from the EMfld rf cavity onto bunch location.
        subroutine getaxfldE_EMfld(z,this,ez1,ezp1,ezpp1)
        implicit none
        include 'mpif.h'
        type (EMfld), intent(in) :: this
        double precision, intent(in) :: z
        double precision, intent(out) :: ez1,ezp1,ezpp1
        double precision:: zz,hstep,slope,zedge,escale
        integer :: klo,khi,k
        integer :: my_rank,ierr

        zedge = this%Param(1)
        escale = this%Param(2)
        zz=z-zedge
        klo=1
        khi=Ndata
1       if(khi-klo.gt.1) then
          k=(khi+klo)/2
          if(zdat(k).gt.zz)then
             khi=k
          else
             klo=k
          endif
          goto 1
        endif
        hstep=zdat(khi)-zdat(klo)
        slope=(edat(khi)-edat(klo))/hstep
        ez1 =edat(klo)+slope*(zz-zdat(klo))
        slope=(epdat(khi)-epdat(klo))/hstep
        ezp1=epdat(klo)+slope*(zz-zdat(klo))
        slope=(eppdat(khi)-eppdat(klo))/hstep
        ezpp1=eppdat(klo)+slope*(zz-zdat(klo))
        ez1 = ez1*escale
        ezp1 = ezp1*escale
        ezpp1 = ezpp1*escale

!        call MPI_COMM_RANK(MPI_COMM_WORLD,my_rank,ierr)
!        if(my_rank.eq.1) then
!         write(19,101)z,ez1/1.0e6,ezp1/1.0e6,zz,float(m),float(init)
!        endif
!        if(my_rank.eq.1) then
!         write(19,101)z,gt
!        endif
  101   format(6(1x,1pe13.6))

        end subroutine getaxfldE_EMfld

!-----------------------------------------------------------------
! get external RF field on axis from analytical function
! Here, we have used a Fouier expansion representation of external field.
! Users should supply the field function as they want.
        subroutine  getaxfldEfc_EMfld(z,this,ez1,ezp1,ezpp1)
        implicit none
        include 'mpif.h'
        double precision, intent(in) :: z
        type (EMfld), intent(in) :: this
        double precision, intent(out) :: ez1,ezp1,ezpp1
        double precision :: zedge,escale,len,zz
        double precision :: pi,zmid
        integer :: i

        pi = 2*asin(1.0)

        zedge = this%Param(1)
        escale = this%Param(2)
        len = this%Length
        zmid = zedge + len/2

        zz=z-zmid
        ez1 = Fcoef(1)/2
        ezp1 = 0.0
        ezpp1 = 0.0
        !Ndata has to be odd number
        do i = 2, (Ndata-1)/2+1
          ez1 = ez1 + Fcoef(2*i-2)*cos((i-1)*2*pi*zz/len) + &
                Fcoef(2*i-1)*sin((i-1)*2*pi*zz/len)
          ezp1 = ezp1 + (i-1)*2*pi/len*(-Fcoef(2*i-2)*sin((i-1)*2*pi*zz/len)+&
                Fcoef(2*i-1)*cos((i-1)*2*pi*zz/len))
          ezpp1 = ezpp1+((i-1)*2*pi/len)**2*(-Fcoef(2*i-2)*cos((i-1)*2*pi*zz/len)&
                  - Fcoef(2*i-1)*sin((i-1)*2*pi*zz/len))
        enddo
        ez1=ez1*escale
        ezp1=ezp1*escale
        ezpp1=ezpp1*escale

        end subroutine getaxfldEfc_EMfld

!-----------------------------------------------------------------
! get external field Ex, Ey, Ez, Bx, Bx, Bz at given position x, y, z, t from
! analytical function. Here we have used Fourier expansion of function. The
! user should supply his own analytical function if needed.
        subroutine  getfld_EMfld(pos,extfld,this)
        implicit none
        include 'mpif.h'
        double precision, dimension(4), intent(in) :: pos
        type (EMfld), intent(in) :: this
        double precision, dimension(6), intent(out) :: extfld
        double precision :: zedge,escale,theta0,ww,len,tt,zz,xl
        double precision :: ez1,ezp1,ezpp1,ezppp,f1,f1p,clite,tmpsin,tmpcos,pi
        double precision :: tmpex,tmpey,tmpez,tmpbx,tmpby,tmpbz,r2
        double precision :: c5,c6,s5,s6,zmid
        integer :: i

        clite = 299792458.e0
        pi = 2*asin(1.0)

        extfld = 0.0
        zedge = this%Param(1)
        escale = this%Param(2)
        len = this%Length
        zmid = zedge + len/2
        ww = this%Param(3)/Scfreq
        xl = Scxl/ww  !real frequency has to be used here
        theta0 = this%Param(4)*asin(1.0)/90
        !move into the tank local coordinate.
        !zz=pos(3)-zedge
!-----------------------------------------------------------------
! get external RF field on axis from analytical function
! Here, we have used a Fouier expansion representation of external field.
! Users should supply the field function as they want.
        zz=pos(3)-zmid
        ez1 = Fcoef(1)/2
        ezp1 = 0.0
        ezpp1 = 0.0
        ezppp = 0.0
        !Ndata has to be odd number
        do i = 2, (Ndata-1)/2+1
          ez1 = ez1 + Fcoef(2*i-2)*cos((i-1)*2*pi*zz/len) + &
                Fcoef(2*i-1)*sin((i-1)*2*pi*zz/len)
          ezp1 = ezp1 + (i-1)*2*pi/len*(-Fcoef(2*i-2)*sin((i-1)*2*pi*zz/len)+&
                Fcoef(2*i-1)*cos((i-1)*2*pi*zz/len))
          ezpp1 = ezpp1+((i-1)*2*pi/len)**2*(-Fcoef(2*i-2)*cos((i-1)*2*pi*zz/len)&
                  - Fcoef(2*i-1)*sin((i-1)*2*pi*zz/len))
          ezppp = ezppp+((i-1)*2*pi/len)**3*(Fcoef(2*i-2)*sin((i-1)*2*pi*zz/len)&
                  - Fcoef(2*i-1)*cos((i-1)*2*pi*zz/len))
        enddo
        ez1=ez1*escale
        ezp1=ezp1*escale
        ezpp1=ezpp1*escale
        ezppp = ezppp*escale
!-----------------------------------------------------------------
! calculate the E and B field using the analytical ez0 on the axis.
! This is based on the assumption of azmuthal symmetry and TM mode.
        tt = pos(4)
        tmpcos = cos(ww*tt+theta0)
        tmpsin = sin(ww*tt+theta0)
        f1 = -(ezpp1+ez1/xl/xl)/4
        f1p = -(ezppp+ezp1/xl/xl)/4
        r2 = pos(1)**2+pos(2)**2
        tmpex = -pos(1)*(ezp1/2+f1p*r2/4)*tmpcos
        tmpey = -pos(2)*(ezp1/2+f1p*r2/4)*tmpcos
        tmpez = (ez1+f1*r2)*tmpcos
        tmpbx = pos(2)/(xl*clite)*(ez1/2+f1*r2/4)*tmpsin
        tmpby = -pos(1)/(xl*clite)*(ez1/2+f1*r2/4)*tmpsin
        tmpbz = 0.0

        extfld(1) = extfld(1) + tmpex
        extfld(2) = extfld(2) + tmpey
        extfld(3) = extfld(3) + tmpez
        extfld(4) = extfld(4) + tmpbx
        extfld(5) = extfld(5) + tmpby
        extfld(6) = extfld(6) + tmpbz

        end subroutine getfld_EMfld

        !get external field with displacement and rotation errors.
        subroutine  getflderrold_EMfld(pos,extfld,this)
        implicit none
        include 'mpif.h'
        double precision, dimension(4), intent(in) :: pos
        type (EMfld), intent(in) :: this
        double precision, dimension(6), intent(out) :: extfld
        double precision :: zedge,escale,theta0,ww,len,tt,zz,xl
        double precision :: ez1,ezp1,ezpp1,ezppp,f1,f1p,clite,tmpsin,tmpcos,pi
        double precision :: tmpex,tmpey,tmpez,tmpbx,tmpby,tmpbz,r2
        double precision :: zmid,dx,dy,anglex,angley,anglez
        double precision, dimension(3) :: temp,tmp
        integer :: i

        clite = 299792458.e0
        pi = 2*asin(1.0)

        zedge = this%Param(1)
        escale = this%Param(2)
        len = this%Length
        zmid = zedge + len/2
        ww = this%Param(3)/Scfreq
        xl = Scxl/ww  !real frequency has to be used here
        theta0 = this%Param(4)*asin(1.0)/90

        dx = this%Param(8)
        dy = this%Param(9)
        anglex = this%Param(10)
        angley = this%Param(11)
        anglez = this%Param(12)

        !move into the tank local coordinate.
        !zz=pos(3)-zedge
        !zz=pos(3)-zmid

        temp(1) = pos(1) - dx
        temp(2) = pos(2) - dy
        tmp(1) = temp(1)*cos(anglez) + temp(2)*sin(anglez)
        tmp(2) = -temp(1)*sin(anglez) + temp(2)*cos(anglez)
        tmp(3) = pos(3)-zmid
        temp(1) = tmp(1)*cos(angley)+tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = -tmp(1)*sin(angley)+tmp(3)*cos(angley)
        tmp(1) = temp(1)
        tmp(2) = temp(2)*cos(anglex)+temp(3)*sin(anglex)
        tmp(3) = -temp(2)*sin(anglex)+temp(3)*cos(anglex)
        zz = tmp(3)

!-----------------------------------------------------------------
! get external RF field on axis from analytical function
! Here, we have used a Fouier expansion representation of external field.
! Users should supply the field function as they want.
        ez1 = Fcoef(1)/2
        ezp1 = 0.0
        ezpp1 = 0.0
        ezppp = 0.0
        !Ndata has to be odd number
        do i = 2, (Ndata-1)/2+1
          ez1 = ez1 + Fcoef(2*i-2)*cos((i-1)*2*pi*zz/len) + &
                Fcoef(2*i-1)*sin((i-1)*2*pi*zz/len)
          ezp1 = ezp1 + (i-1)*2*pi/len*(-Fcoef(2*i-2)*sin((i-1)*2*pi*zz/len)+&
                Fcoef(2*i-1)*cos((i-1)*2*pi*zz/len))
          ezpp1 = ezpp1+((i-1)*2*pi/len)**2*(-Fcoef(2*i-2)*cos((i-1)*2*pi*zz/len)&
                  - Fcoef(2*i-1)*sin((i-1)*2*pi*zz/len))
          ezppp = ezppp+((i-1)*2*pi/len)**3*(Fcoef(2*i-2)*sin((i-1)*2*pi*zz/len)&
                  - Fcoef(2*i-1)*cos((i-1)*2*pi*zz/len))
        enddo
        ez1=ez1*escale
        ezp1=ezp1*escale
        ezpp1=ezpp1*escale
        ezppp = ezppp*escale
!-----------------------------------------------------------------
! calculate the E and B field using the analytical ez0 on the axis.
! This is based on the assumption of azmuthal symmetry and TM mode.
        tt = pos(4)
        tmpcos = cos(ww*tt+theta0)
        tmpsin = sin(ww*tt+theta0)
        f1 = -(ezpp1+ez1/xl/xl)/4
        f1p = -(ezppp+ezp1/xl/xl)/4
        r2 = tmp(1)**2+tmp(2)**2
        extfld(1) = -tmp(1)*(ezp1/2+f1p*r2/4)*tmpcos
        extfld(2) = -tmp(2)*(ezp1/2+f1p*r2/4)*tmpcos
        extfld(3) = (ez1+f1*r2)*tmpcos
        extfld(4) = tmp(2)/(xl*clite)*(ez1/2+f1*r2/4)*tmpsin
        extfld(5) = -tmp(1)/(xl*clite)*(ez1/2+f1*r2/4)*tmpsin
        extfld(6) = 0.0

        !for E field
        tmp(1) = extfld(1)
        tmp(2) = extfld(2)*cos(anglex)-extfld(3)*sin(anglex)
        tmp(3) = extfld(2)*sin(anglex)+extfld(3)*cos(anglex)
        temp(1) = tmp(1)*cos(angley)-tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = tmp(1)*sin(angley)+tmp(3)*cos(angley)
        extfld(1) = temp(1)*cos(anglez) - temp(2)*sin(anglez)
        extfld(2) = temp(1)*sin(anglez) + temp(2)*cos(anglez)
        extfld(3) = temp(3)

        !for B field
        tmp(1) = extfld(4)
        tmp(2) = extfld(5)*cos(anglex)-extfld(6)*sin(anglex)
        tmp(3) = extfld(5)*sin(anglex)+extfld(6)*cos(anglex)
        temp(1) = tmp(1)*cos(angley)-tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = tmp(1)*sin(angley)+tmp(3)*cos(angley)
        extfld(4) = temp(1)*cos(anglez) - temp(2)*sin(anglez)
        extfld(5) = temp(1)*sin(anglez) + temp(2)*cos(anglez)
        extfld(6) = temp(3)

        end subroutine getflderrold_EMfld

        !get external field with displacement and rotation errors.
        subroutine  getflderr_EMfld(pos,extfld,this,dx,dy,anglex,angley,anglez)
        implicit none
        include 'mpif.h'
        double precision, dimension(4), intent(in) :: pos
        type (EMfld), intent(in) :: this
        double precision, intent(in) :: dx,dy,anglex,angley,anglez
        double precision, dimension(6), intent(out) :: extfld
        double precision :: zedge,escale,theta0,ww,len,tt,zz,xl
        double precision :: ez1,ezp1,ezpp1,ezppp,f1,f1p,clite,tmpsin,tmpcos,pi
        double precision :: tmpex,tmpey,tmpez,tmpbx,tmpby,tmpbz,r2
        double precision :: zmid
        double precision, dimension(3) :: temp,tmp
        integer :: i

        clite = 299792458.e0
        pi = 2*asin(1.0)

        zedge = this%Param(1)
        escale = this%Param(2)
        len = this%Length
        zmid = zedge + len/2
        ww = this%Param(3)/Scfreq
        xl = Scxl/ww  !real frequency has to be used here
        theta0 = this%Param(4)*asin(1.0)/90

        !move into the tank local coordinate.
        !zz=pos(3)-zedge
        !zz=pos(3)-zmid

        temp(1) = pos(1) - dx
        temp(2) = pos(2) - dy
        tmp(1) = temp(1)*cos(anglez) + temp(2)*sin(anglez)
        tmp(2) = -temp(1)*sin(anglez) + temp(2)*cos(anglez)
        tmp(3) = pos(3) - zmid
        temp(1) = tmp(1)*cos(angley)+tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = -tmp(1)*sin(angley)+tmp(3)*cos(angley)
        tmp(1) = temp(1)
        tmp(2) = temp(2)*cos(anglex)+temp(3)*sin(anglex)
        tmp(3) = -temp(2)*sin(anglex)+temp(3)*cos(anglex)
        zz = tmp(3)

!-----------------------------------------------------------------
! get external RF field on axis from analytical function
! Here, we have used a Fouier expansion representation of external field.
! Users should supply the field function as they want.
        ez1 = Fcoef(1)/2
        ezp1 = 0.0
        ezpp1 = 0.0
        ezppp = 0.0
        !Ndata has to be odd number
        do i = 2, (Ndata-1)/2+1
          ez1 = ez1 + Fcoef(2*i-2)*cos((i-1)*2*pi*zz/len) + &
                Fcoef(2*i-1)*sin((i-1)*2*pi*zz/len)
          ezp1 = ezp1 + (i-1)*2*pi/len*(-Fcoef(2*i-2)*sin((i-1)*2*pi*zz/len)+&
                Fcoef(2*i-1)*cos((i-1)*2*pi*zz/len))
          ezpp1 = ezpp1+((i-1)*2*pi/len)**2*(-Fcoef(2*i-2)*cos((i-1)*2*pi*zz/len)&
                  - Fcoef(2*i-1)*sin((i-1)*2*pi*zz/len))
          ezppp = ezppp+((i-1)*2*pi/len)**3*(Fcoef(2*i-2)*sin((i-1)*2*pi*zz/len)&
                  - Fcoef(2*i-1)*cos((i-1)*2*pi*zz/len))
        enddo
        ez1=ez1*escale
        ezp1=ezp1*escale
        ezpp1=ezpp1*escale
        ezppp = ezppp*escale
!-----------------------------------------------------------------
! calculate the E and B field using the analytical ez0 on the axis.
! This is based on the assumption of azmuthal symmetry and TM mode.
        tt = pos(4)
        tmpcos = cos(ww*tt+theta0)
        tmpsin = sin(ww*tt+theta0)
        f1 = -(ezpp1+ez1/xl/xl)/4
        f1p = -(ezppp+ezp1/xl/xl)/4
        r2 = tmp(1)**2+tmp(2)**2
        extfld(1) = -tmp(1)*(ezp1/2+f1p*r2/4)*tmpcos
        extfld(2) = -tmp(2)*(ezp1/2+f1p*r2/4)*tmpcos
        extfld(3) = (ez1+f1*r2)*tmpcos
        extfld(4) = tmp(2)/(xl*clite)*(ez1/2+f1*r2/4)*tmpsin
        extfld(5) = -tmp(1)/(xl*clite)*(ez1/2+f1*r2/4)*tmpsin
        extfld(6) = 0.0

        !for E field
        tmp(1) = extfld(1)
        tmp(2) = extfld(2)*cos(anglex)-extfld(3)*sin(anglex)
        tmp(3) = extfld(2)*sin(anglex)+extfld(3)*cos(anglex)
        temp(1) = tmp(1)*cos(angley)-tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = tmp(1)*sin(angley)+tmp(3)*cos(angley)
        extfld(1) = temp(1)*cos(anglez) - temp(2)*sin(anglez)
        extfld(2) = temp(1)*sin(anglez) + temp(2)*cos(anglez)
        extfld(3) = temp(3)

        !for B field
        tmp(1) = extfld(4)
        tmp(2) = extfld(5)*cos(anglex)-extfld(6)*sin(anglex)
        tmp(3) = extfld(5)*sin(anglex)+extfld(6)*cos(anglex)
        temp(1) = tmp(1)*cos(angley)-tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = tmp(1)*sin(angley)+tmp(3)*cos(angley)
        extfld(4) = temp(1)*cos(anglez) - temp(2)*sin(anglez)
        extfld(5) = temp(1)*sin(anglez) + temp(2)*cos(anglez)
        extfld(6) = temp(3)

        end subroutine getflderr_EMfld

!-----------------------------------------------------------------
! get the discrete Er,Etheta,Ez, Br, Btheta, Bz as a function or
! "r" at given "z".
        subroutine getfld6_EMfld(this,z,extfld6)
        type (EMfld), intent(in) :: this
        double precision, intent(in) :: z
        double precision, dimension(6,NrIntvRf+1), intent(out) :: extfld6
        double precision :: zz,hz,ef

        hz = (ZmaxRf-ZminRf)/NzIntvRf
        zz = z-this%Param(1)
        iz = zz/hz + 1
        if(iz.eq.(NzIntvRf+1)) then
           iz = iz - 1
        endif
        iz1 = iz+1
        ef = (iz*hz - zz)/hz

! get the external field along "r" at given "z" by linear interpolation
! Here, we have only Er, Ez, and Btheta for TM mode.
        do j = 1, NrIntvRf+1
          extfld6(1,j) = erdata(iz,j)*ef + erdata(iz1,j)*(1.0-ef)
          extfld6(2,j) = 0.0
          extfld6(3,j) = ezdata(iz,j)*ef + ezdata(iz1,j)*(1.0-ef)
          extfld6(4,j) = 0.0
          extfld6(5,j) = btdata(iz,j)*ef + btdata(iz1,j)*(1.0-ef)
          extfld6(6,j) = 0.0
        enddo

        end subroutine getfld6_EMfld

        subroutine getfld6err_EMfld(this,pos,extfld,dx,dy,anglex,angley,&
                                    anglez)
        type (EMfld), intent(in) :: this
        double precision, dimension(4), intent(in) :: pos
        double precision, intent(in) :: dx,dy,anglex,angley,anglez
        double precision, dimension(6), intent(out) :: extfld
        double precision :: zz,hz,efz,hr,efr,len,zedge,escale,ww,theta0,&
                            tt,er,ez,bt,tmpcos,tmpsin
        double precision, dimension(3) :: tmp,temp
        integer :: iz,iz1,ir,ir1

        len = this%Length
        zedge = this%Param(1)
        escale = this%Param(2)
        ww = this%Param(3)/Scfreq
        theta0 = this%Param(4)*asin(1.0)/90

        hz = (ZmaxRf-ZminRf)/NzIntvRf
        hr = (RmaxRf-RminRf)/NrIntvRf

        temp(1) = pos(1) - dx
        temp(2) = pos(2) - dy
        tmp(1) = temp(1)*cos(anglez) + temp(2)*sin(anglez)
        tmp(2) = -temp(1)*sin(anglez) + temp(2)*cos(anglez)
        tmp(3) = pos(3) - zedge
        temp(1) = tmp(1)*cos(angley)+tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = -tmp(1)*sin(angley)+tmp(3)*cos(angley)
        tmp(1) = temp(1)
        tmp(2) = temp(2)*cos(anglex)+temp(3)*sin(anglex)
        tmp(3) = -temp(2)*sin(anglex)+temp(3)*cos(anglex)
        zz = tmp(3)
        tt = pos(4)
        tmpcos = cos(ww*tt+theta0)*escale
        tmpsin = sin(ww*tt+theta0)*escale

        if(zz.le.0) then
          iz = 1
          efz = zz/hz
        else if(zz.ge.(ZmaxRf-ZminRf)) then
          iz = NzIntvRf
          efz = (zz-(iz-1)*hz)/hz
        else
          iz = zz/hz + 1
          efz = (zz-(iz-1)*hz)/hz
        endif
        iz1 = iz+1

        rr = sqrt(tmp(1)*tmp(1)+tmp(2)*tmp(2))
        ir = rr/hr + 1
        if(rr.eq.0.0) then
          rr = 1.0e-10
        endif
        if(ir.ge.NrIntvRf) then
           print*,"ir: ",ir,rr,tmp(1),tmp(2)
           stop
        endif
        efr = (rr-(ir-1)*hr)/hr
        ir1 = ir+1

        ! get the external field along "r" at given "z" by linear interpolation
        ! Here, we have only Er, Ez, and Btheta for TM mode.
        er = (erdata(iz,ir)*(1.0-efz)*(1.0-efr)+erdata(iz1,ir)*efz*(1.0-efr)+&
             erdata(iz1,ir1)*efz*efr+erdata(iz,ir1)*(1.0-efz)*efr)*tmpcos
        ez = (ezdata(iz,ir)*(1.0-efz)*(1.0-efr)+ezdata(iz1,ir)*efz*(1.0-efr)+&
             ezdata(iz1,ir1)*efz*efr+ezdata(iz,ir1)*(1.0-efz)*efr)*tmpcos
        bt = (btdata(iz,ir)*(1.0-efz)*(1.0-efr)+btdata(iz1,ir)*efz*(1.0-efr)+&
             btdata(iz1,ir1)*efz*efr+btdata(iz,ir1)*(1.0-efz)*efr)*tmpsin

        extfld(1) = er*tmp(1)/rr
        extfld(2) = er*tmp(2)/rr
        extfld(3) = ez
        extfld(4) = -bt*tmp(2)/rr
        extfld(5) = bt*tmp(1)/rr
        extfld(6) = 0.0

        !for E field
        tmp(1) = extfld(1)
        tmp(2) = extfld(2)*cos(anglex)-extfld(3)*sin(anglex)
        tmp(3) = extfld(2)*sin(anglex)+extfld(3)*cos(anglex)
        temp(1) = tmp(1)*cos(angley)-tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = tmp(1)*sin(angley)+tmp(3)*cos(angley)
        extfld(1) = temp(1)*cos(anglez) - temp(2)*sin(anglez)
        extfld(2) = temp(1)*sin(anglez) + temp(2)*cos(anglez)
        extfld(3) = temp(3)

        !for B field
        tmp(1) = extfld(4)
        tmp(2) = extfld(5)*cos(anglex)-extfld(6)*sin(anglex)
        tmp(3) = extfld(5)*sin(anglex)+extfld(6)*cos(anglex)
        temp(1) = tmp(1)*cos(angley)-tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = tmp(1)*sin(angley)+tmp(3)*cos(angley)
        extfld(4) = temp(1)*cos(anglez) - temp(2)*sin(anglez)
        extfld(5) = temp(1)*sin(anglez) + temp(2)*cos(anglez)
        extfld(6) = temp(3)

        end subroutine getfld6err_EMfld

!-----------------------------------------------------------------
! get the discrete Ex,Ey,Ez, Bx, By, Bz as a function of x and y at
! given "z".
        subroutine getfld6xyz_EMfld(this,z,extfld6xyz)
        type (EMfld), intent(in) :: this
        double precision, intent(in) :: z
        double precision, dimension(6,NxIntvRfg+1,NyIntvRfg+1),&
                          intent(out) :: extfld6xyz
        double precision :: zz,hz,ef

        hz = (ZmaxRfg-ZminRfg)/NzIntvRfg
        zz = z-this%Param(1)
        iz = zz/hz + 1
        if(iz.eq.(NzIntvRfg+1)) then
           iz = iz - 1
        endif
        iz1 = iz+1
        ef = (iz*hz - zz)/hz

! get the external field along "r" at given "z" by linear interpolation
! Here, we have only Er, Ez, and Btheta for TM mode.
        do j = 1, NyIntvRfg+1
          do i = 1, NxIntvRfg+1
            extfld6xyz(1,i,j) = Exgrid(i,j,iz)*ef + Exgrid(i,j,iz1)*(1.0-ef)
            extfld6xyz(2,i,j) = Eygrid(i,j,iz)*ef + Eygrid(i,j,iz1)*(1.0-ef)
            extfld6xyz(3,i,j) = Ezgrid(i,j,iz)*ef + Ezgrid(i,j,iz1)*(1.0-ef)
            extfld6xyz(4,i,j) = Bxgrid(i,j,iz)*ef + Bxgrid(i,j,iz1)*(1.0-ef)
            extfld6xyz(5,i,j) = Bygrid(i,j,iz)*ef + Bygrid(i,j,iz1)*(1.0-ef)
            extfld6xyz(6,i,j) = Bzgrid(i,j,iz)*ef + Bzgrid(i,j,iz1)*(1.0-ef)
          enddo
        enddo

        end subroutine getfld6xyz_EMfld

        subroutine getfld6xyzerr_EMfld(this,pos,extfld,dx,dy,anglex,angley,&
                                    anglez)
        type (EMfld), intent(in) :: this
        double precision, dimension(4), intent(in) :: pos
        double precision, intent(in) :: dx,dy,anglex,angley,anglez
        double precision, dimension(6), intent(out) :: extfld
        double precision :: zz,hz,efz,hx,efx,hy,efy,len,zedge,escale,ww,theta0,&
                        tt,exn,eyn,ezn,bxn,byn,bzn,tmpcos,tmpsin
        double precision, dimension(3) :: tmp,temp
        integer :: iz,iz1,ix,ix1,iy,iy1

        len = this%Length
        zedge = this%Param(1)
        escale = this%Param(2)
        ww = this%Param(3)/Scfreq
        theta0 = this%Param(4)*asin(1.0)/90

        hx = (XmaxRfg-XminRfg)/NxIntvRfg
        hy = (YmaxRfg-YminRfg)/NyIntvRfg
        hz = (ZmaxRfg-ZminRfg)/NzIntvRfg

        temp(1) = pos(1) - dx
        temp(2) = pos(2) - dy
        tmp(1) = temp(1)*cos(anglez) + temp(2)*sin(anglez)
        tmp(2) = -temp(1)*sin(anglez) + temp(2)*cos(anglez)
        tmp(3) = pos(3) - zedge
        temp(1) = tmp(1)*cos(angley)+tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = -tmp(1)*sin(angley)+tmp(3)*cos(angley)
        tmp(1) = temp(1)
        tmp(2) = temp(2)*cos(anglex)+temp(3)*sin(anglex)
        tmp(3) = -temp(2)*sin(anglex)+temp(3)*cos(anglex)
        zz = tmp(3)
        tt = pos(4)
        tmpcos = cos(ww*tt+theta0)*escale
        tmpsin = sin(ww*tt+theta0)*escale

        if(tmp(1).le.XminRfg) then
          ix = 1
          efx = (XminRfg-tmp(1)+ix*hx)/hx
        else if(tmp(1).ge.XmaxRfg) then
          ix = NxIntvRfg
          efx = (XminRfg-tmp(1)+ix*hx)/hx
        else
          ix = (tmp(1)-XminRfg)/hx + 1
          efx = (XminRfg-tmp(1)+ix*hx)/hx
        endif
        ix1 = ix+1

        if(tmp(2).le.YminRfg) then
          iy = 1
          efy = (YminRfg-tmp(2)+iy*hy)/hy
        else if(tmp(2).ge.YmaxRfg) then
          iy = NyIntvRfg
          efy = (YminRfg-tmp(2)+iy*hy)/hy
        else
          iy = (tmp(2)-YminRfg)/hy + 1
          efy = (YminRfg-tmp(2)+iy*hy)/hy
        endif
        iy1 = iy+1

        if(zz.le.0.0) then
          iz = 1
          efz = 1.0 - zz/hz
        else if(zz.ge.(ZmaxRfg-ZminRfg)) then
          iz = NzIntvRfg
          efz = 1.0 - (zz-(iz-1)*hz)/hz
        else
          iz = zz/hz + 1
          efz = 1.0 - (zz-(iz-1)*hz)/hz
        endif
        iz1 = iz+1

        exn = Exgrid(ix,iy,iz)*efx*efy*efz + &
          Exgrid(ix,iy1,iz)*efx*(1-efy)*efz  + &
          Exgrid(ix,iy1,iz1)*efx*(1-efy)*(1-efz) + &
          Exgrid(ix,iy,iz1)*efx*efy*(1-efz) + &
          Exgrid(ix1,iy,iz1)*(1-efx)*efy*(1-efz) + &
          Exgrid(ix1,iy1,iz1)*(1-efx)*(1-efy)*(1-efz) + &
          Exgrid(ix1,iy1,iz)*(1-efx)*(1-efy)*efz + &
          Exgrid(ix1,iy,iz)*(1-efx)*efy*efz
        eyn = Eygrid(ix,iy,iz)*efx*efy*efz + &
          Eygrid(ix,iy1,iz)*efx*(1-efy)*efz  + &
          Eygrid(ix,iy1,iz1)*efx*(1-efy)*(1-efz) + &
          Eygrid(ix,iy,iz1)*efx*efy*(1-efz) + &
          Eygrid(ix1,iy,iz1)*(1-efx)*efy*(1-efz) + &
          Eygrid(ix1,iy1,iz1)*(1-efx)*(1-efy)*(1-efz) + &
          Eygrid(ix1,iy1,iz)*(1-efx)*(1-efy)*efz + &
          Eygrid(ix1,iy,iz)*(1-efx)*efy*efz
        ezn = Ezgrid(ix,iy,iz)*efx*efy*efz + &
          Ezgrid(ix,iy1,iz)*efx*(1-efy)*efz  + &
          Ezgrid(ix,iy1,iz1)*efx*(1-efy)*(1-efz) + &
          Ezgrid(ix,iy,iz1)*efx*efy*(1-efz) + &
          Ezgrid(ix1,iy,iz1)*(1-efx)*efy*(1-efz) + &
          Ezgrid(ix1,iy1,iz1)*(1-efx)*(1-efy)*(1-efz) + &
          Ezgrid(ix1,iy1,iz)*(1-efx)*(1-efy)*efz + &
          Ezgrid(ix1,iy,iz)*(1-efx)*efy*efz
        bxn = Bxgrid(ix,iy,iz)*efx*efy*efz + &
          Bxgrid(ix,iy1,iz)*efx*(1-efy)*efz  + &
          Bxgrid(ix,iy1,iz1)*efx*(1-efy)*(1-efz) + &
          Bxgrid(ix,iy,iz1)*efx*efy*(1-efz) + &
          Bxgrid(ix1,iy,iz1)*(1-efx)*efy*(1-efz) + &
          Bxgrid(ix1,iy1,iz1)*(1-efx)*(1-efy)*(1-efz) + &
          Bxgrid(ix1,iy1,iz)*(1-efx)*(1-efy)*efz + &
          Bxgrid(ix1,iy,iz)*(1-efx)*efy*efz
        byn = Bygrid(ix,iy,iz)*efx*efy*efz + &
          Bygrid(ix,iy1,iz)*efx*(1-efy)*efz  + &
          Bygrid(ix,iy1,iz1)*efx*(1-efy)*(1-efz) + &
          Bygrid(ix,iy,iz1)*efx*efy*(1-efz) + &
          Bygrid(ix1,iy,iz1)*(1-efx)*efy*(1-efz) + &
          Bygrid(ix1,iy1,iz1)*(1-efx)*(1-efy)*(1-efz) + &
          Bygrid(ix1,iy1,iz)*(1-efx)*(1-efy)*efz + &
          Bygrid(ix1,iy,iz)*(1-efx)*efy*efz
        bzn = Bzgrid(ix,iy,iz)*efx*efy*efz + &
          Bzgrid(ix,iy1,iz)*efx*(1-efy)*efz  + &
          Bzgrid(ix,iy1,iz1)*efx*(1-efy)*(1-efz) + &
          Bzgrid(ix,iy,iz1)*efx*efy*(1-efz) + &
          Bzgrid(ix1,iy,iz1)*(1-efx)*efy*(1-efz) + &
          Bzgrid(ix1,iy1,iz1)*(1-efx)*(1-efy)*(1-efz) + &
          Bzgrid(ix1,iy1,iz)*(1-efx)*(1-efy)*efz + &
          Bzgrid(ix1,iy,iz)*(1-efx)*efy*efz

        extfld(1) = exn*tmpcos
        extfld(2) = eyn*tmpcos
        extfld(3) = ezn*tmpcos
        extfld(4) = bxn*tmpsin
        extfld(5) = byn*tmpsin
        extfld(6) = bzn*tmpsin

        !transform
        !for E field
        tmp(1) = extfld(1)
        tmp(2) = extfld(2)*cos(anglex)-extfld(3)*sin(anglex)
        tmp(3) = extfld(2)*sin(anglex)+extfld(3)*cos(anglex)
        temp(1) = tmp(1)*cos(angley)-tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = tmp(1)*sin(angley)+tmp(3)*cos(angley)
        extfld(1) = temp(1)*cos(anglez) - temp(2)*sin(anglez)
        extfld(2) = temp(1)*sin(anglez) + temp(2)*cos(anglez)
        extfld(3) = temp(3)

        !for B field
        tmp(1) = extfld(4)
        tmp(2) = extfld(5)*cos(anglex)-extfld(6)*sin(anglex)
        tmp(3) = extfld(5)*sin(anglex)+extfld(6)*cos(anglex)
        temp(1) = tmp(1)*cos(angley)-tmp(3)*sin(angley)
        temp(2) = tmp(2)
        temp(3) = tmp(1)*sin(angley)+tmp(3)*cos(angley)
        extfld(4) = temp(1)*cos(anglez) - temp(2)*sin(anglez)
        extfld(5) = temp(1)*sin(anglez) + temp(2)*cos(anglez)
        extfld(6) = temp(3)

        end subroutine getfld6xyzerr_EMfld

      end module EMfldclass

max_stars_repo_name: hhslepicka/IMPACT-Z

max_issues_repo_name: hhslepicka/IMPACT-Z

max_forks_repo_name: hhslepicka/IMPACT-Z


text:
*> \brief \b SROUNDUP_LWORK
*
*  =========== DOCUMENTATION ===========
*
* Online html documentation available at
*            http://www.netlib.org/lapack/explore-html/
*
*  Definition:
*  ===========
*
*      REAL             FUNCTION SROUNDUP_LWORK( LWORK )
*
*     .. Scalar Arguments ..
*      INTEGER          LWORK
*     ..
*
*
*> \par Purpose:
*  =============
*>
*> \verbatim
*>
*> SROUNDUP_LWORK deals with a subtle bug with returning LWORK as a Float.
*> This routine guarantees it is rounded up instead of down by
*> multiplying LWORK by 1+eps when it is necessary, where eps is the relative machine precision.
*> E.g.,
*>
*>        float( 16777217            ) == 16777216
*>        float( 16777217 ) * (1.+eps) == 16777218
*>
*> \return SROUNDUP_LWORK
*> \verbatim
*>         SROUNDUP_LWORK >= LWORK.
*>         SROUNDUP_LWORK is guaranteed to have zero decimal part.
*> \endverbatim
*
*  Arguments:
*  ==========
*
*> \param[in] LWORK Workspace size.
*
*  Authors:
*  ========
*
*> \author Weslley Pereira, University of Colorado Denver, USA
*
*> \ingroup auxOTHERauxiliary
*
*> \par Further Details:
*  =====================
*>
*> \verbatim
*>  This routine was inspired in the method `magma_zmake_lwork` from MAGMA.
*>  \see https://bitbucket.org/icl/magma/src/master/control/magma_zauxiliary.cpp
*> \endverbatim
*
*  =====================================================================
      REAL             FUNCTION SROUNDUP_LWORK( LWORK )
*
*  -- LAPACK auxiliary routine --
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
*     .. Scalar Arguments ..
      INTEGER          LWORK
*     ..
*
* =====================================================================
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC         EPSILON, REAL, INT
*     ..
*     .. Executable Statements ..
*     ..
      SROUNDUP_LWORK = REAL( LWORK )
*
      IF( INT( SROUNDUP_LWORK ) .LT. LWORK ) THEN
*         Force round up of LWORK
          SROUNDUP_LWORK = SROUNDUP_LWORK * ( 1.0E+0 + EPSILON(0.0E+0) )
      ENDIF
*
      RETURN
*
*     End of SROUNDUP_LWORK
*
      END

max_stars_repo_name: quellyn/lapack

max_issues_repo_name: quellyn/lapack

max_forks_repo_name: quellyn/lapack