!
! This file is part of SACAMOS, State of the Art CAble MOdels in Spice. 
! It was developed by the University of Nottingham and the Netherlands Aerospace 
! Centre (NLR) for ESA under contract number 4000112765/14/NL/HK.
! 
! Copyright (C) 2016-2017 University of Nottingham
! 
! SACAMOS is free software: you can redistribute it and/or modify it under the 
! terms of the GNU General Public License as published by the Free Software 
! Foundation, either version 3 of the License, or (at your option) any later 
! version.
! 
! SACAMOS is distributed in the hope that it will be useful, but 
! WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 
! or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License 
! for more details.
! 
! A copy of the GNU General Public License version 3 can be found in the 
! file GNU_GPL_v3 in the root or at <http://www.gnu.org/licenses/>.
! 
! SACAMOS uses the EISPACK library (in /SRC/EISPACK). EISPACK is subject to 
! the GNU Lesser General Public License. A copy of the GNU Lesser General Public 
! License version can be found in the file GNU_LGPL in the root of EISPACK 
! (/SRC/EISPACK ) or at <http://www.gnu.org/licenses/>.
! 
! The University of Nottingham can be contacted at: ggiemr@nottingham.ac.uk
!
! File Contents:
! SUBROUTINE gt_zero_check(p,cable_spec_error,cable_name,message)
! SUBROUTINE cylindrical_check(r,cable_spec_error,cable_name,message)
! SUBROUTINE cylindrical_with_dielectric_check(r,rd,cable_spec_error,cable_name,message)
! SUBROUTINE coax_with_dielectric_check(r,rs,rd,cable_spec_error,cable_name,message)
! SUBROUTINE twisted_pair_check(r,rd,s,cable_spec_error,cable_name,message)
! SUBROUTINE shielded_twisted_pair_check(r,rd,s,rs,rd2,cable_spec_error,cable_name,message)
! SUBROUTINE spacewire_check(r,rd,s,rs,rd2,stpr,rs2,rd3,cable_spec_error,cable_name,message)
! SUBROUTINE rectangular_check(w,h,cable_spec_error,cable_name,message)
! SUBROUTINE flex_cable_check(w,h,s,dxo,dyo,cable_spec_error,cable_name,message)
! SUBROUTINE dielectric_check(epsr,cable_spec_error,cable_name,message)
! SUBROUTINE transfer_impedance_check(Zt,cable_spec_error,cable_name,message)
! SUBROUTINE conductivity_check(sigma,cable_spec_error,cable_name,message)
! SUBROUTINE surface_impedance_check(Zt,sigma,rs,t,cable_spec_error,cable_name,message)
!
! NAME
!     gt_zero_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that a parameter,p, is greater than zero
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 14/11/2016 CJS 
!     3/8/2017 return with error message

!     16/11/2017 CJS Include network synthesis process to replace s-domain transfer functions

!
!
SUBROUTINE gt_zero_check(p,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  real(dp),intent(IN) :: p  ! input test parameter
  
  logical,intent(INOUT)  :: cable_spec_error    ! error flag
  
  character(LEN=line_length),intent(IN)              :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

if (p.LE.0d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Parameter is less than or equal to zero'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

RETURN

END SUBROUTINE gt_zero_check
!
! NAME
!     cylindrical_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining a cylindrical cable are consistent
!     The wire radius must be greater than zero
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!
SUBROUTINE cylindrical_check(r,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  real(dp),intent(IN) :: r
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

if (r.LE.0d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Radius is less than or equal to zero'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

RETURN

END SUBROUTINE cylindrical_check
!
! NAME
!     cylindrical_with_dielectric_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining a cylindrical cable with dielectric are consistent
!     the wire radius must be greater than zero and 
!     the dielectric radius must be greater than or equal to the conductor radius.
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!
SUBROUTINE cylindrical_with_dielectric_check(r,rd,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  real(dp),intent(IN) :: r
  real(dp),intent(IN) :: rd
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

CALL cylindrical_check(r,cable_spec_error,cable_name,message)

if (rd.LT.r) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Dielectric radius is less than the condcutor radius'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

RETURN

END SUBROUTINE cylindrical_with_dielectric_check
!
! NAME
!     coax_with_dielectric_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining a coax cable are consistent
!     the inner wire radius must be greater than zero and 
!     the shield radius must be greater than the conductor radius.
!     the dielectric radius must be greater than or equal to the shield radius.
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!     
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!
SUBROUTINE coax_with_dielectric_check(r,rs,rd,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  real(dp),intent(IN) :: r
  real(dp),intent(IN) :: rs
  real(dp),intent(IN) :: rd
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

CALL cylindrical_check(r,cable_spec_error,cable_name,message)

if (rs.LE.r) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Shield radius is less than or equal to the condcutor radius'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

if (rd.LT.rs) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Dielectric radius is less than the shield radius'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

RETURN

END SUBROUTINE coax_with_dielectric_check
!
! NAME
!     twisted_pair_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining a twisted pair cable with dielectric are consistent
!     the wire with dielectric checks must be satisfied plus
!     the two wires must not touch i.e. the twisted pair separation must be greater than twice the dielectric radius
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!
SUBROUTINE twisted_pair_check(r,rd,s,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  real(dp),intent(IN) :: r
  real(dp),intent(IN) :: rd
  real(dp),intent(IN) :: s
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

CALL cylindrical_with_dielectric_check(r,rd,cable_spec_error,cable_name,message)

if (s.LE.rd*2d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='twisted pair separation is less than or equal to twice the dielectric radius'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

RETURN

END SUBROUTINE twisted_pair_check
!
! NAME
!     shielded_twisted_pair_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining a shielded twisted pair cable with dielectric are consistent
!     the twisted pair checks must be satisfied
!     the shield must not intersect the twisted pair wires
!     the shield dielectric radius must be greater than or equal to the shield radiius
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!
SUBROUTINE shielded_twisted_pair_check(r,rd,s,rs,rd2,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  real(dp),intent(IN) :: r
  real(dp),intent(IN) :: rd
  real(dp),intent(IN) :: s
  real(dp),intent(IN) :: rs
  real(dp),intent(IN) :: rd2
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

CALL twisted_pair_check(r,rd,s,cable_spec_error,cable_name,message)

if (rs.LE.s/2d0+rd) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Shield intersects the twisted pair wires'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

if (rd2.LT.rs) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Dielectric radius is less than the shield radius'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

RETURN

END SUBROUTINE shielded_twisted_pair_check
!
! NAME
!     spacewire_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining spacewire are consistent
!     the shielded twisted pair checks must be satisfied
!     the four shielded twisted pairs must not intersect
!     the outer shield must not intersect the shielded twisted pairs
!     the outer dielectric radius must be greater than or equal to the outer shield radius
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!
SUBROUTINE spacewire_check(r,rd,s,rs,rd2,stpr,rs2,rd3,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  real(dp),intent(IN) :: r
  real(dp),intent(IN) :: rd
  real(dp),intent(IN) :: s
  real(dp),intent(IN) :: rs
  real(dp),intent(IN) :: rd2
  real(dp),intent(IN) :: stpr
  real(dp),intent(IN) :: rs2
  real(dp),intent(IN) :: rd3
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

CALL shielded_twisted_pair_check(r,rd,s,rs,rd2,cable_spec_error,cable_name,message)

if (stpr.LE.sqrt(2d0)*rd2) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Shielded twisted pairs intersect'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

if (rs2.LE.rd2+stpr) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Outer shield intersects the shielded twisted pairs'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

if (rd3.LT.rs2) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Outer dielectric radius is less than the outer shield radius'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

RETURN

END SUBROUTINE spacewire_check
!
! NAME
!     rectangular_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining rectangular conductors are consistent
!     the conductor width and height must both be greater than zero
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!
SUBROUTINE rectangular_check(w,h,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  real(dp),intent(IN) :: w
  real(dp),intent(IN) :: h
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

if (w.LE.0d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Conductor width is less than or equal to zero'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

if (h.LE.0d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Conductor height is less than or equal to zero'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

RETURN

END SUBROUTINE rectangular_check
!
! NAME
!     flex_cable_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining flex cable are consistent
!     the conductor width and height must both be greater than zero
!     the conductor separation must be greater than zero if there is more more than one conductor
!     the dielectric offset in x and y must be greater than or equal to zero
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!
SUBROUTINE flex_cable_check(nc,w,h,s,dxo,dyo,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  integer,intent(IN)  :: nc
  real(dp),intent(IN) :: w
  real(dp),intent(IN) :: h
  real(dp),intent(IN) :: s
  real(dp),intent(IN) :: dxo
  real(dp),intent(IN) :: dyo
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

CALL rectangular_check(w,h,cable_spec_error,cable_name,message)

if (s.LT.0d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='conductor separation is less than zero'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

if ((nc.GT.1).AND.(s.EQ.0d0)) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='conductor separation is zero'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

if (dxo.LT.0d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Dielectric offset in x is less than zero'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

if (dyo.LT.0d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Dielectric offset in y is less than zero'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

RETURN

END SUBROUTINE flex_cable_check
!
! NAME
!     dielectric_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining dielectric models are consistent
!     d.c. relative permittivity must be greater than or equal to 1
!     high frequency relative permittivity must be greater than or equal to 1
!     the poles of the rational function must be on LHS of the s-plane (time domain stability requirement)
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!

SUBROUTINE dielectric_check(epsr,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  type(Sfilter),intent(IN) :: epsr
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

real(dp) :: epsr_value
logical  :: stable

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

! Check that the d.c. relative permittivity is greater than or equal to 1.0

epsr_value=epsr%a%coeff(0)/epsr%b%coeff(0)

if (epsr_value.LT.1d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Relative permittivity is less than 1.0 at d.c.'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

! Check that the relative permittivity is greater than or equal to 1.0 at high frequency

epsr_value=evaluate_Sfilter_high_frequency_limit(epsr)

if (epsr_value.LT.1d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Relative permittivity is less than 1.0 at high frequency'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

! check that the poles of the filter function are in the LHS of the s plane

CALL test_filter_pole_stability(epsr,stable)
if (.NOT.stable) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Relative permittivity model has unstable pole(s)'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if

END SUBROUTINE dielectric_check
!
! NAME
!     not_FD_dielectric_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the dielectric model is not frequency dependent as this is a limitation for certain models
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!

SUBROUTINE not_FD_dielectric_check(epsr,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  type(Sfilter),intent(IN) :: epsr
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

logical  :: stable

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

! Check that the d.c. relative permittivity is greater than or equal to 1.0

if ( (epsr%a%order.NE.0).OR.(epsr%b%order.NE.0) ) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Relative permittivity model is frequency dependent'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if


END SUBROUTINE not_FD_dielectric_check
!
! NAME
!     transfer impedance_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining transfer impedance are consistent
!     the d.c. transfer impedance must be greater than or equal to zero
!     the poles of the rational function must be on LHS of the s-plane (time domain stability requirement)     
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!

SUBROUTINE transfer_impedance_check(Zt,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  type(Sfilter),intent(IN) :: ZT

  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

real(dp) :: Rdc
logical  :: stable

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

! Check that the d.c. resistance is greater than or equal to zero

Rdc=ZT%a%coeff(0)/ZT%b%coeff(0)

if (Rdc.LT.0d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Transfer impedance is less than zero at d.c.'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

! check that the poles of the filter function are in the LHS of the s plane

CALL test_filter_pole_stability(ZT,stable)
if (.NOT.stable) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Transfer impedance model has unstable pole(s)'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if

END SUBROUTINE transfer_impedance_check
!
! NAME
!     conductivity_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the conductivity specified is consistent
!     the conductivity specified must be greater than or equal to zero. 
!     Note that zero conductivity indicates a perfect conductor in this software.
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!

SUBROUTINE conductivity_check(sigma,cable_spec_error,cable_name,message)

USE type_specifications

IMPLICIT NONE

! variables passed to subroutine

  real(dp),intent(IN) :: sigma
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

if (sigma.LT.0d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Conductivity is less than zero'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

END SUBROUTINE conductivity_check
!
! NAME
!     surface_impedance_check
!
! AUTHORS
!     Chris Smartt
!
! DESCRIPTION
!     check that the parameters defining surface impedance for cylindrical shields are consistent
!     The thickness,t , must be greater than or equal to zero on input
!     If t=0 then calculate an equivalent thickness such that the surface impedance is equal to the
!     transfer impedance at low frequency i.e. make the surface and transfer impedance models consistent
!     If t=0 and the conductivity=0 and the transfer impedance at d.c. is not zero then the shield model
!     cannot be specified and an error is flagged.
!     If the check fails the cable_spec_error flag is set on return, otherwise it is left unchanged
!
! COMMENTS
!      
!
! HISTORY
!
!     started 28/10/2016 CJS 
!     3/8/2017 return with error message
!
!

SUBROUTINE surface_impedance_check(ZT,sigma,rs,t,cable_spec_error,cable_name,message)

USE type_specifications
USE constants

IMPLICIT NONE

! variables passed to subroutine

  type(Sfilter) :: ZT

  real(dp),intent(INOUT) :: sigma
  real(dp),intent(IN) :: rs
  real(dp),intent(INOUT) :: t
  
  logical,intent(INOUT)  :: cable_spec_error
  
  character(LEN=line_length),intent(IN) :: cable_name  
  character(LEN=error_message_length),intent(INOUT) :: message 

! local variables

  real(dp) :: Rdc

! START

if (cable_spec_error) RETURN ! return if an error has already been flagged

CALL conductivity_check(sigma,cable_spec_error,cable_name,message)

if (t.LT.0d0) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Shield thickness is less than zero'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if 

! calculate Rdc, the transfer impedance at d.c.
 Rdc=ZT%a%coeff(0)/ZT%b%coeff(0)

! If t=0 then calculate an equivalent thickness such that the surface impedance is equal to the
! transfer impedance at low frequency
  
if ((t.EQ.0d0).AND.(sigma.NE.0d0)) then
! we need to calculate the thickness to be consistent with the transfer impedance at d.c. i.e. R_dc = ZT_dc
  
  if (Rdc.NE.0d0) then
! calculate the 'equivalent thickness' for the shield
    t=1d0/(2d0*pi*rs*sigma*Rdc)
  else
! The transfer impedance is zero so set the surface impedance to be zero by setting sigma=0
    t=0d0
    sigma=0d0
  end if
  
end if
  
if ((t.EQ.0d0).AND.(sigma.EQ.0d0).AND.(Rdc.GT.small)) then
  write(*,*)'Error in cable:',trim(cable_name)
  message='Shield conductivity and thickness cannot both be set to zero if there is a non-zero transfer impedance'
  write(*,'(A)')trim(message)
  cable_spec_error=.TRUE.
  RETURN
end if

END SUBROUTINE surface_impedance_check