C
C  User subroutine for defining body current density

      subroutine udecurrent(bodycurrent,predef,coords,nBlock,
     $     i_array,niarray,r_array,nrarray,c_array,ncarray)
c
      include 'aba_param.inc'
c
      dimension bodycurrent(nBlock,*),i_array(*),r_array(*)
      dimension predef(nBlock,2,*),coords(nBlock,*)
      character*80 c_array(*)
c
      parameter(i_udecurr_kstep = 1,
     $          i_udecurr_kinc  = 2,
     $          i_udecurr_noel  = 3,
     $          i_udecurr_npt   = 4,
     $          i_udecurr_jltyp = 5,
     $          i_udecurr_phase = 6,
     $          i_udecurr_proc  = 7,
     $          i_udecurr_nfld  = 8)
c
      parameter(ir_udecurr_time_1   = 1,
     $          ir_udecurr_time_2   = 2)
c
      parameter(i_jltyp_cj = 1,
     $          i_jltyp_r  = 2)
c
      parameter(i_proc_lf_th = 1)
c
      parameter(i_udecurr_phase_real = 1,
     $          i_udecurr_phase_imag = 2)
c
      parameter (zero    = 0.d0,
     $           one     = 1.d0,
     $           sigma   = 0.58d0,
     $           pi      = 3.141592653589793d0,
     $           two     = 2.d0,
     $           sixty   = 60.d0)
c
c     r = radius
c     E_phi = unit vector in circumferential direction
c     i is the standard notation used for imaginary numbers
c     A = r E_phi * (1 + i)
c     J = i*omega*sigma*A
c
      coord1 = coords(1,1)
      coord2 = coords(1,2)
      radius = sqrt(coord1*coord1 + coord2*coord2)
c
      bodycurrent(1,1)=zero
      bodycurrent(1,3)=zero
c
      if (i_array(i_udecurr_phase).eq.i_udecurr_phase_real) then
        tmp = -two*pi*sixty*sigma
      else if (i_array(i_udecurr_phase).eq.i_udecurr_phase_imag) then
        tmp = two*pi*sixty*sigma
      end if
c
      bodycurrent(1,2)=tmp*radius
c
      return
      end