! 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 .
!
! 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 .
!
! The University of Nottingham can be contacted at: ggiemr@nottingham.ac.uk
!
! File Contents:
! SUBROUTINE write_incident_field_circuit
!
! NAME
! write_incident_field_circuit
!
! AUTHORS
! Chris Smartt
!
! DESCRIPTION
! This code writes the circuit components required to implement the incident field excitation model
! The Spice model is described in Theory_Manual_Section 3.8
! and the Spice circuit is seen in Theory_Manual_Figures 3.17, and 3.18
!
! INPUTS REQUIRED
! 1. nodes for connection of the incident field sources
! 2. time delays for the viictim domain modes
! 3. bundle length
! 4. modal decomposition matrices
! 5. victim domain propagation correction filter functions
! 6. pre-calculated coefficients related to the incident plane wave field excitation
! 7. nodes for the specification of the incident field function
! 8. next free node and reference node for the spice sub-circuit
!
!
! OUTPUTS
! The components required to implement the model are written to
! the subcircuit file
!
! COMMENTS
! Write the circuit elements required for the incident field excitation model
!
! The excitation function comes from the voltage between two subcircuit termination nodes.
!
! HISTORY
!
! STAGE 6 developments started 16/06/2016 CJS.
! Include ground plane 28/6/2016 CJS
! Replace Tz transmission lines with a very small delay by a very small resistance - causes long runtimes in Pspice. 6/12/2016.
! 9/5/2017 CJS Document software with reference to Theory_Manual
!
SUBROUTINE write_incident_field_circuit(n_victim_domain_modes, &
Vv_end1_node,Vv_end2_node, &
T_victim, &
length, &
TVI_victim,TI_victim, &
Hpv_filter, &
Tz,alpha0,alphaL, &
next_free_node,Einc_node1,Einc_node2,vref_node,victim_domain )
USE type_specifications
USE general_module
USE constants
USE spice_cable_bundle_module
USE filter_module
USE maths
IMPLICIT NONE
! variables passed to the subroutine
integer,intent(IN) :: n_victim_domain_modes ! number of modes in the victim domain
! victim domain node lists for connection of the incident field source terms at each end.
integer,intent(IN) :: Vv_end1_node(1:n_victim_domain_modes)
integer,intent(IN) :: Vv_end2_node(1:n_victim_domain_modes)
! mode delays for the victim domain transmission line modes
real(dp),intent(IN) :: T_victim(1:n_victim_domain_modes)
! modal decomposition matrices for the victim domain
real(dp),intent(IN) :: TVI_victim(1:n_victim_domain_modes,1:n_victim_domain_modes)
real(dp),intent(IN) :: TI_victim(1:n_victim_domain_modes,1:n_victim_domain_modes)
real(dp),intent(IN) :: length ! bundle length
TYPE(Sfilter),intent(IN) :: Hpv_filter(1:n_victim_domain_modes) ! victim mode proopagation correction terms
real(dp),intent(IN) :: Tz ! incident field excitation delay time
real(dp),intent(IN) :: alpha0(1:n_victim_domain_modes) ! pre-calculated model constants
real(dp),intent(IN) :: alphaL(1:n_victim_domain_modes) ! pre-calculated model constants
integer,intent(INOUT) :: next_free_node ! the next free node number for the spice sub-circuit
integer,intent(IN) :: Einc_node1 ! excitation function sub-circuit terminal node 1
integer,intent(IN) :: Einc_node2 ! excitation function sub-circuit terminal node 2
integer,intent(IN) :: vref_node ! spice sub-circuit reference node number
integer,intent(IN) :: victim_domain ! victim domain number
! local variables
integer :: v_mode
real(dp) :: Rlarge=1D10 ! large resistance
real(dp) :: Z0_delay=50d0 ! transmission line impedance for pure delay lines
! Names for delay lines
character(len=spice_name_length) :: delay_line_Tz_name
character(len=spice_name_length) :: delay_line_Tv_name(1:n_victim_domain_modes)
character(len=spice_name_length) :: delay_line_TzPTv_name(1:n_victim_domain_modes)
character(len=spice_name_length) :: delay_line_ZC_Tz_name
character(len=spice_name_length) :: delay_line_ZC_Tv_name(1:n_victim_domain_modes)
character(len=spice_name_length) :: delay_line_ZC_TzPTv_name(1:n_victim_domain_modes)
character(len=spice_name_length) :: delay_line_E1_Tz_name
character(len=spice_name_length) :: delay_line_E1_Tv_name(1:n_victim_domain_modes)
character(len=spice_name_length) :: delay_line_E1_TzPTv_name(1:n_victim_domain_modes)
character(len=spice_name_length) :: E_Einc_l_name(1:n_victim_domain_modes)
character(len=spice_name_length) :: E_Einc_s_name(1:n_victim_domain_modes)
character(len=spice_name_length) :: E_Einc_filter_name(1:n_victim_domain_modes)
! names for special case Tz=Tvictim components
character(len=spice_name_length) :: G_Einc_derivative_name(1:n_victim_domain_modes)
character(len=spice_name_length) :: L_Einc_derivative_name(1:n_victim_domain_modes)
! working strings
character(len=spice_name_length) :: name1,name2,name3
character(len=spice_name_length) :: Einc_string
! String used if a transmission line delay is very small and the T element is replaced by a resistance
character(len=spice_name_length) :: Rtempstring
! circuit to combine incident field terms
character(len=spice_name_length) :: combine_delays_s_E_name(1:n_victim_domain_modes,1:2)
character(len=spice_name_length) :: R_combine_delays_s_name(1:n_victim_domain_modes)
character(len=spice_name_length) :: combine_delays_l_E_name(1:n_victim_domain_modes,1:2)
character(len=spice_name_length) :: R_combine_delays_l_name(1:n_victim_domain_modes)
real(dp) :: Evalue
! Nodes for delay lines
! delay line Tz
integer :: delay_line_Tz_s_nodes
integer :: delay_line_Tz_l_nodes
! delay line Tv
integer :: delay_line_Tv_s_nodes(1:n_victim_domain_modes)
integer :: delay_line_Tv_l_nodes(1:n_victim_domain_modes)
! delay line Tz+Tv
integer :: delay_line_TzPTv_s_nodes(1:n_victim_domain_modes)
integer :: delay_line_TzPTv_l_nodes(1:n_victim_domain_modes)
integer :: last_Hjw_source_s_node(1:n_victim_domain_modes)
integer :: last_Hjw_source_l_node(1:n_victim_domain_modes)
integer :: combine_delays_s_Enode
integer :: combine_delays_l_Enode
real(dp) :: gain ! gain value for s-domain transfer functions
real(dp) :: Tz_minus_Tv ! difference between source and victim mode velocities
logical :: Tz_equal_Tv ! if incident field delay and victim mode delay are the same then we have a different circuit topology
! loop variables
integer :: row,i
! nodes for special case Tz=Tvictim components
integer :: Einc_derivative_node(1:n_victim_domain_modes)
! START
if (verbose) write(*,*)'CALLED write_incient_field_excitation_circuit'
! check for the case where Tz<0 and stop with an error for now.
! Eventually we will cope with this situation by reversing the model.
if (Tz.LT.0d0) then
run_status='ERROR in write_incient_field_excitation_circuit. Tz is less than zero'
CALL write_program_status()
STOP 1
end if
CALL write_spice_comment('START OF INCIDENT FIELD EXCITATION MODELS')
! Set Tz delay line node numbers
! delay line nodes for positive z propagation, source end
next_free_node=next_free_node+1
delay_line_Tz_s_nodes=next_free_node
! delay line nodes for positive z propagation, load end
next_free_node=next_free_node+1
delay_line_Tz_l_nodes=next_free_node
! Set Tz delay line component names
delay_line_Tz_name='T_Tz_Einc'
! mode impedance
delay_line_ZC_Tz_name='RZC_Tz_Einc'
! source terms
delay_line_E1_Tz_name='E1_Tz_Einc'
! Write Tz delay line components
! Theory_Manual_Fig 3.17, Theory_Manual_Eqn 3.168b
if (Tz.GT.Tz_min_delay) then
! significant delay so use a T element
CALL write_spice_comment('Incident field delay lines, Tz delay')
write(spice_model_file_unit,'(A30,4I6,A4,E16.6,A4,E16.6)')delay_line_Tz_name, &
delay_line_Tz_s_nodes,vref_node, &
delay_line_Tz_l_nodes,vref_node, &
' Z0=',Z0_delay,' TD=',Tz
else
! the delay is very small so just use a very small series resistance instead of a delay line
CALL write_spice_comment('Incident field delay lines, Tz delay=0. Delay line replaced by small resistance')
Rtempstring='R'//trim(delay_line_Tz_name)
write(spice_model_file_unit,'(A30,2I6,E16.6)')Rtempstring, &
delay_line_Tz_s_nodes,delay_line_Tz_l_nodes,Rsmall
end if
! modal impedances on modal delay lines, Tz
CALL write_spice_comment('Matched impedance: Tz delay')
write(spice_model_file_unit,'(A30,2I6,E16.6)')delay_line_ZC_Tz_name, &
delay_line_Tz_l_nodes,vref_node,Z0_delay
! delay line controlled source for positive z propagation, Ts
CALL write_spice_comment('Delay line controlled sources Tz delay')
write(spice_model_file_unit,'(A30,4I6,E16.6)')delay_line_E1_Tz_name,&
delay_line_Tz_s_nodes,vref_node, &
Einc_node1,Einc_node2,1.0
! End of Tz delay line
do v_mode=1,n_victim_domain_modes ! loop over victim domain modes
! create Einc_string which labels the transfer impedance model number plus the source mode and victim mode numbers
name1='EINC_vm_'
CALL add_integer_to_string(name1,v_mode,Einc_string)
! we always need the Tz+T_victim(v_mode) delay lines so write these components now
! Theory_Manual_Fig 3.17, Theory_Manual_Eqn 3.168a
! Set TzPTv delay line nodes
! delay line nodes, source end
next_free_node=next_free_node+1
delay_line_TzPTv_s_nodes(v_mode)=next_free_node
! delay line nodes, load end
next_free_node=next_free_node+1
delay_line_TzPTv_l_nodes(v_mode)=next_free_node
! Set TzPTv delay line component names
delay_line_TzPTv_name(v_mode)='T_TzPTv_'//trim(Einc_string)
! mode impedance
delay_line_ZC_TzPTv_name(v_mode)='RZC_TzPTv_'//trim(Einc_string)
! source terms
delay_line_E1_TzPTv_name(v_mode)='E1_TzPTv_'//trim(Einc_string)
! Write TzPTv delay lines
! Theory_Manual_Fig 3.17, Theory_Manual_Eqn 3.168a
CALL write_spice_comment('Delay line Tz+T_victim(v_mode)')
write(spice_model_file_unit,'(A30,4I6,A4,E16.6,A4,E16.6)')delay_line_TzPTv_name(v_mode),&
delay_line_TzPTv_s_nodes(v_mode),vref_node, &
delay_line_TzPTv_l_nodes(v_mode),vref_node, &
' Z0=',Z0_delay,' TD=',Tz+T_victim(v_mode)
! modal impedances on modal delay lines, TzPTv
CALL write_spice_comment('Matched impedance: Tz+T_victim(v_mode) delay')
write(spice_model_file_unit,'(A30,2I6,E16.6)')delay_line_ZC_TzPTv_name(v_mode), &
delay_line_TzPTv_l_nodes(v_mode),vref_node,Z0_delay
! delay line controlled source, TzPTv
CALL write_spice_comment('Incident field delay line controlled source: Tz+T_victim(v_mode) delay')
write(spice_model_file_unit,'(A30,4I6,E16.6)')delay_line_E1_TzPTv_name(v_mode),&
delay_line_TzPTv_s_nodes(v_mode),vref_node, &
Einc_node1,Einc_node2,1.0
! End of TzPTv delay lines
! Test for the special case when the source and victim mode delays are the same (or very close)
! In this case we need to use a different model to avoid a singularity in the normal model
! Theory_Manual_Section 3.8.1, Theory_Manual_Eqn 3.169
Tz_minus_Tv=Tz-T_victim(v_mode)
if (abs(Tz_minus_Tv).GT.Einc_min_delay) then
Tz_equal_Tv=.FALSE.
! The whole incident field excitation circuit is implemented using delay lines
! Theory_Manual_Fig 3.17
! Set Tv delay line nodes
! delay line nodes source end
next_free_node=next_free_node+1
delay_line_Tv_s_nodes(v_mode)=next_free_node
! delay line nodes load end
next_free_node=next_free_node+1
delay_line_Tv_l_nodes(v_mode)=next_free_node
! Set Tv delay line component names
delay_line_Tv_name(v_mode)='T_Tv_'//trim(Einc_string)
! mode impedance
delay_line_ZC_Tv_name(v_mode)='RZC_Tv_'//trim(Einc_string)
! source terms
delay_line_E1_Tv_name(v_mode)='E1_Tv_'//trim(Einc_string)
! Write Tv delay lines
CALL write_spice_comment('Incident field delay line T_victim(v_mode)')
write(spice_model_file_unit,'(A30,4I6,A4,E16.6,A4,E16.6)')delay_line_Tv_name(v_mode), &
delay_line_Tv_s_nodes(v_mode),vref_node, &
delay_line_Tv_l_nodes(v_mode),vref_node, &
' Z0=',Z0_delay,' TD=',T_victim(v_mode)
! modal impedances on modal delay lines, T_victim(v_mode)
CALL write_spice_comment('Matched impedance T_victim(v_mode) delay')
write(spice_model_file_unit,'(A30,2I6,E16.6)')delay_line_ZC_Tv_name(v_mode), &
delay_line_Tv_l_nodes(v_mode),vref_node,Z0_delay
! delay line controlled source for positive z propagation, T_victim(v_mode)
CALL write_spice_comment('Controlled source for T_victim(v_mode) delay')
write(spice_model_file_unit,'(A30,4I6,E16.6)')delay_line_E1_Tv_name(v_mode) &
,delay_line_Tv_s_nodes(v_mode),vref_node, &
Einc_node1,Einc_node2,1.0
else
! Special case Tz=T_victim(v_mode)as see in Theory_Manual_Eqn 3.169 Theory_Manual_Fig 3.18
Tz_equal_Tv=.TRUE.
! **** The special case required time derivative circuits operating on the delayed incident field excitation****
! new node for time derivative delayed incident field function
next_free_node=next_free_node+1
Einc_derivative_node(v_mode)=next_free_node
! ****** Names for the special case circuit for Tsource=Tvictim
G_Einc_derivative_name(v_mode)='G_Vp_ddt_'//trim(Einc_string)
L_Einc_derivative_name(v_mode)='L_Vp_ddt_'//trim(Einc_string)
! ***** Inductive circuit to calculate the time derivative of delayed Einc *****
CALL write_spice_comment('Controlled source for derivative of delayed Einc i.e. Einc(0,t-Tz)')
write(spice_model_file_unit,'(A30,4I6,E16.6)')G_Einc_derivative_name(v_mode), &
Einc_derivative_node(v_mode),vref_node, &
delay_line_Tz_l_nodes,vref_node,1.0
CALL write_spice_comment('1H inductor for derivative of delayed Einc i.e. Einc(0,t-Tz)')
write(spice_model_file_unit,'(A30,2I6,E16.6)')L_Einc_derivative_name(v_mode), &
Einc_derivative_node(v_mode),vref_node,1.0
end if ! Special case Tz-Tv =0
! The remaining part of the circuit combines all the contributions to the victim mode voltage source
! So as to evaluate the incident field terms in Theory_Manual_Eqn 3.168a,b
! We create the nodes for the summation circuit as we go
! create Einc_string which labels the transfer impedance model number plus victim mode number
name1='EINC_vm_'
CALL add_integer_to_string(name1,v_mode,Einc_string)
! first add the contributions which are common to both forms of circuit
! Incident field excitation voltage source names
combine_delays_s_E_name(v_mode,1)='E_zt_dsum_s_'//trim(Einc_string)//'_E1'
combine_delays_s_E_name(v_mode,2)='E_zt_dsum_s_'//trim(Einc_string)//'_E2'
! START OF CIRCUIT TO COMBINE INCIDENT FIELD EXCITATION TERMS
! calculation of V_victim at z=0. Theory_Manual_Eqn 3.168a
CALL write_spice_comment('Circuit to combine incident field excitation terms')
! Einc, no delay: Theory_Manual_Eqn 3.168a, line 2, term 1.
next_free_node=next_free_node+1
Evalue=+alpha0(v_mode)/(Tz+T_victim(v_mode))
write(spice_model_file_unit,'(A30,4I6,E16.6)')combine_delays_s_E_name(v_mode,1) &
,next_free_node,vref_node &
,Einc_node1,Einc_node2 &
,Evalue
combine_delays_s_Enode=next_free_node
! Einc, delay=T_victim+Tz: Theory_Manual_Eqn 3.168a, line 2, term 2
next_free_node=next_free_node+1
Evalue=-alpha0(v_mode)/(Tz+T_victim(v_mode))
write(spice_model_file_unit,'(A30,4I6,E16.6)')combine_delays_s_E_name(v_mode,2) &
,next_free_node,combine_delays_s_Enode &
,delay_line_TzPTv_l_nodes(v_mode),vref_node &
,Evalue
combine_delays_s_Enode=next_free_node
! calculation of V_victim at z=L. Theory_Manual_Eqn 3.168 a
if (.NOT.Tz_equal_Tv) then
! normal form based on delay lines, Theory_Manual_Eqn 3.168b
combine_delays_l_E_name(v_mode,1)='E_zt_dsum_l_'//trim(Einc_string)//'_E1'
combine_delays_l_E_name(v_mode,2)='E_zt_dsum_l_'//trim(Einc_string)//'_E2'
! Einc, delay=T_victim: Theory_Manual_Eqn 3.168b, line 2, term 1.
next_free_node=next_free_node+1
Evalue=+alphaL(v_mode)/(T_victim(v_mode)-Tz)
write(spice_model_file_unit,'(A30,4I6,E16.6)')combine_delays_l_E_name(v_mode,1) &
,next_free_node,vref_node &
,delay_line_Tv_l_nodes(v_mode),vref_node &
,Evalue
combine_delays_l_Enode=next_free_node
! Einc, delay=Tz: Theory_Manual_Eqn 3.168b, line 2, term 2.
next_free_node=next_free_node+1
Evalue=-alphaL(v_mode)/(T_victim(v_mode)-Tz)
write(spice_model_file_unit,'(A30,4I6,E16.6)')combine_delays_l_E_name(v_mode,2) &
,next_free_node,combine_delays_l_Enode &
,delay_line_Tz_l_nodes,vref_node &
,Evalue
combine_delays_l_Enode=next_free_node
else
! we need to calcuate the contribution to z=L using the time derivative of the delayed incident field
! Theory_Manual_Eqn 3.169, Theory_Manual_Fig 3.18
combine_delays_l_E_name(v_mode,1)='E_zt_dsum_l_'//trim(Einc_string)//'_E1'
next_free_node=next_free_node+1
Evalue=-alphaL(v_mode)
write(spice_model_file_unit,'(A30,4I6,E16.6)')combine_delays_l_E_name(v_mode,1) &
,next_free_node,combine_delays_l_Enode &
,Einc_derivative_node(v_mode),vref_node &
,Evalue
combine_delays_l_Enode=next_free_node
end if ! special case T_victim=Tz
CALL write_spice_comment('Large resistance to complete the circuit for the series voltage sources')
! Large resistance to complete the circuit for the series voltage sources
R_combine_delays_l_name(v_mode)='R_Einc_dsum_l_'//trim(Einc_string)
write(spice_model_file_unit,'(A30,2I6,E16.6)')R_combine_delays_l_name(v_mode),combine_delays_l_Enode,vref_node,Rlarge
! Large resistance to complete the circuit for the series voltage sources
R_combine_delays_s_name(v_mode)='R_Einc_dsum_s_'//trim(Einc_string)
write(spice_model_file_unit,'(A30,2I6,E16.6)')R_combine_delays_s_name(v_mode),combine_delays_s_Enode,vref_node,Rlarge
! Write filter function for propagation correction filter
! Theory_Manual_Equation 3.170
CALL write_spice_comment('Incident field excitation sources, end 1')
E_Einc_s_name(v_mode)='Einc_s_'//trim(Einc_string)
CALL write_s_domain_controlled_voltage_source(E_Einc_s_name(v_mode), &
combine_delays_s_Enode,vref_node, &
Vv_end1_node(v_mode),Vref_node, &
Hpv_filter(v_mode),1d0,vref_node,next_free_node) ! note: gain set to 1.0
CALL write_spice_comment('Incident field excitation sources, end 2')
E_Einc_l_name(v_mode)='Einc_l_'//trim(Einc_string)
CALL write_s_domain_controlled_voltage_source(E_Einc_l_name(v_mode), &
combine_delays_l_Enode,vref_node, &
Vv_end2_node(v_mode),Vref_node, &
Hpv_filter(v_mode),1d0,vref_node,next_free_node) ! note: gain set to 1.0
end do ! next victim mode
RETURN
END SUBROUTINE write_incident_field_circuit