! ! This file is part of SACAMOS, State of the Art CAble MOdels for 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-2018 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: ! MODULE cable_bundle_module ! CONTAINS ! SUBROUTINE read_cable_bundle ! SUBROUTINE write_cable_bundle ! SUBROUTINE deallocate_cable_bundle ! SUBROUTINE check_cables_wrt_ground_plane ! SUBROUTINE check_cable_intersections ! ! NAME ! MODULE cable_bundle_module ! ! Data strunctures and subroutines relating to cable bundles ! ! COMMENTS ! Need to improve error handling in file reading ! Not sure that some arrays in the type_specification need to be in this structure as they are only used in one subroutine... ! ! HISTORY ! started 25/11/2015 CJS ! stage 3 developments started 21/3/2015 CJS ! generalise conductor impedance for transfer impedance model (stage 5) CJS 12/5/2016 ! CJS 25/08/2016 updates for revised transfer impedance and conductor impedance model for shields ! 16/11/2017 CJS Include network synthesis process to replace s-domain transfer functions ! MODULE cable_bundle_module USE type_specifications USE cable_module IMPLICIT NONE ! This is a general structure used to hold conductor based lists, TYPE::conductor_list_type integer :: n_elements integer,allocatable :: element(:) END TYPE conductor_list_type ! This is the main structure used to hold all the bundle information TYPE::bundle_specification_type character(LEN=line_length) :: version character(LEN=line_length) :: bundle_name integer :: tot_n_conductors integer :: tot_n_shielded_conductors integer :: tot_n_external_conductors ! dimension of the global matrix system (total number of conductors-1) integer :: system_dimension ! cable information, total number of cables including and excluding the ground plane (if it exists) integer :: n_cables integer :: n_cables_without_ground_plane ! Array of cable specification data for all the cables in this bundle type(cable_specification_type),allocatable :: cable(:) ! list of cable positions in the bundle cross section and rotation angles real(dp),allocatable :: cable_x_offset(:) real(dp),allocatable :: cable_y_offset(:) real(dp),allocatable :: cable_angle(:) ! ground plane information: note we have two notations at the moment. ! This will need to be sorted out and one only adopted logical :: ground_plane_present real(dp) :: ground_plane_angle real(dp) :: ground_plane_offset real(dp) :: ground_plane_sigma ! ground plane conductivity for FastHenry2 real(dp) :: ground_plane_w ! ground plane width for FastHenry2 real(dp) :: ground_plane_h ! ground plane height for FastHenry2 real(dp) :: ground_plane_Rdc ! ground plane dc resistance for FastHenry2 integer :: ground_plane_nsegx ! ground plane number of segments in x for FastHenry2 integer :: ground_plane_nsegz ! ground plane number of segments in z for FastHenry2 integer :: ground_plane_nh ! ground plane number of layers in h for FastHenry2 real(dp) :: ground_plane_x,ground_plane_y,ground_plane_theta real(dp) :: ground_plane_nx,ground_plane_ny integer :: ground_plane_cable_side ! +1 if the cables are on the normal direction side of the ground plane, -1 otherwise ! domain based information integer :: tot_n_domains integer,allocatable :: n_conductors(:) ! number of conductors in each domain type(conductor_list_type),allocatable :: terminal_conductor_list(:) ! external conductor number for coonductors in each domain type(matrix),allocatable :: L(:) ! 'high frequency' inductance matrix for each domain type(matrix),allocatable :: C(:) ! 'high frequency' capacitance matrix for each domain type(Sfilter_matrix),allocatable :: Z(:) ! frequency dependent impedance matrix for each domain type(Sfilter_matrix),allocatable :: Y(:) ! frequency dependent admittance matrix for each domain ! global conductor based information type(matrix) :: global_MI ! global domain current transformation matrix type(matrix) :: global_MV ! global domain voltage transformation matrix type(matrix) :: global_L ! global 'high frequency' inductance matrix type(matrix) :: global_C ! global 'high frequency' capacitance matrix type(Sfilter_matrix) :: global_Z ! global frequency dependent impedance matrix type(Sfilter_matrix) :: global_Y ! global frequency dependent admittance matrix ! global conductor impedance model information for each conductor type(conductor_impedance_model),allocatable :: conductor_impedance(:) ! global conductor position in the bundle cross-section real(dp),allocatable :: conductor_x_offset(:) real(dp),allocatable :: conductor_y_offset(:) ! arrays for cross referecing numbering information from the global (external) conductor numbering and the domain based numbering systems logical,allocatable :: terminal_conductor_is_shield_flag(:) integer,allocatable :: terminal_conductor_to_inner_domain(:) integer,allocatable :: terminal_conductor_to_outer_domain(:) integer,allocatable :: terminal_conductor_to_global_domain_conductor(:) integer,allocatable :: terminal_conductor_to_local_domain_conductor(:) integer,allocatable :: terminal_conductor_to_reference_terminal_conductor(:) ! Y matrix element function fitting information for frequency dependent dielectrics ! model order for filter fitting integer :: Y_fit_model_order ! frequency range specification type(frequency_specification) :: Y_fit_freq_spec character(LEN=line_length),allocatable :: conductor_label(:) END TYPE bundle_specification_type CONTAINS ! The following subroutines apply to all bundles and cable types ! NAME ! SUBROUTINE read_cable_bundle(unit) ! ! read cable bundle structure from a specified unit ! ! COMMENTS ! ! ! HISTORY ! started 2/12/2015 CJS ! stage 3 developments started 21/3/2015 CJS ! generalise conductor impedance for transfer impedance model (stage 5) CJS 12/5/2016 ! SUBROUTINE read_cable_bundle(bundle,file_unit) USE type_specifications USE constants USE general_module USE cable_module USE maths IMPLICIT NONE ! variables passed to subroutine type(bundle_specification_type),intent(INOUT) :: bundle integer ,intent(IN) :: file_unit ! local variables character(len=filename_length) :: filename logical :: file_exists character(len=line_length) :: line integer :: cable integer :: n_domains integer :: domain integer :: matrix_dimension integer :: n_conductors integer :: conductor integer :: ierr ! START ! Open the (.bundle) file to read filename=trim(MOD_bundle_lib_dir)//trim(bundle%bundle_name)//bundle_file_extn inquire(file=trim(filename),exist=file_exists) if (.NOT.file_exists) then run_status='ERROR, Cannot find the required bundle file:'//trim(filename) CALL write_program_status() STOP 1 end if open(unit=file_unit,file=trim(filename)) if (verbose) write(*,*)'Opened file:',trim(filename) ! read .bundle file read(file_unit,'(A)',ERR=9000)bundle%version read(file_unit,'(A)',ERR=9000)bundle%bundle_name ! read cable information read(file_unit,*,ERR=9000)bundle%n_cables_without_ground_plane read(file_unit,*,ERR=9000)bundle%n_cables ! allocate and read cable information and cable positions in bundle ALLOCATE( bundle%cable(1:bundle%n_cables) ) ALLOCATE( bundle%cable_x_offset(1:bundle%n_cables) ) ALLOCATE( bundle%cable_y_offset(1:bundle%n_cables) ) ALLOCATE( bundle%cable_angle(1:bundle%n_cables) ) do cable=1,bundle%n_cables_without_ground_plane read(file_unit,'(A)')bundle%cable(cable)%cable_name CALL read_cable(bundle%cable(cable),cable_file_unit) read(file_unit,*,ERR=9000)bundle%cable_x_offset(cable),bundle%cable_y_offset(cable),bundle%cable_angle(cable) ! convert angle to radians bundle%cable_angle(cable)=bundle%cable_angle(cable)*pi/180d0 end do ! next cable ! read ground plane specification read(file_unit,'(A)',IOSTAT=ierr)line if (ierr.NE.0) then run_status='ERROR reading ground plane present/ absent information' CALL write_program_status() STOP 1 end if CALL convert_to_lower_case(line,line_length) bundle%ground_plane_present=(line.eq.'ground_plane') if (bundle%ground_plane_present) then read(file_unit,*,IOSTAT=ierr)bundle%ground_plane_x,bundle%ground_plane_y,bundle%ground_plane_angle read(file_unit,*,IOSTAT=ierr)bundle%ground_plane_nx,bundle%ground_plane_ny read(file_unit,*,IOSTAT=ierr)bundle%ground_plane_cable_side if (ierr.NE.0) then run_status='ERROR reading ground plane position and angle' CALL write_program_status() STOP 1 end if ! convert ground plane angle to radians bundle%ground_plane_angle=bundle%ground_plane_angle*pi/180d0 ! calculate offset bundle%ground_plane_offset=-bundle%ground_plane_x*sin(bundle%ground_plane_angle)+ & bundle%ground_plane_y*cos(bundle%ground_plane_angle) read(file_unit,*,IOSTAT=ierr)bundle%ground_plane_sigma read(file_unit,*,IOSTAT=ierr)bundle%ground_plane_w read(file_unit,*,IOSTAT=ierr)bundle%ground_plane_h read(file_unit,*,IOSTAT=ierr)bundle%ground_plane_Rdc read(file_unit,*,IOSTAT=ierr)bundle%ground_plane_nsegx read(file_unit,*,IOSTAT=ierr)bundle%ground_plane_nsegz read(file_unit,*,IOSTAT=ierr)bundle%ground_plane_nh if (ierr.NE.0) then run_status='ERROR reading ground plane sigma, w h and Rdc' CALL write_program_status() STOP 1 end if ! set the geometric data for the last cable in the list i.e. the ground plane cable=bundle%n_cables bundle%cable_x_offset(cable)=bundle%ground_plane_x bundle%cable_y_offset(cable)=bundle%ground_plane_y bundle%cable_angle(cable)=bundle%ground_plane_angle-pi/2d0 ! note there are two ground plane spec formats at the moment... else ! set x, y, angle and offset to 0 bundle%ground_plane_x=0d0 bundle%ground_plane_y=0d0 bundle%ground_plane_nx=0d0 bundle%ground_plane_ny=0d0 bundle%ground_plane_cable_side=0 bundle%ground_plane_angle=0d0 bundle%ground_plane_offset=0d0 end if ! read system dimension information read(file_unit,*)bundle%tot_n_conductors read(file_unit,*)bundle%tot_n_external_conductors read(file_unit,*)bundle%system_dimension ! read domain information read(file_unit,*,ERR=9000)bundle%tot_n_domains ALLOCATE( bundle%n_conductors(1:bundle%tot_n_domains) ) ALLOCATE( bundle%L(1:bundle%tot_n_domains) ) ALLOCATE( bundle%C(1:bundle%tot_n_domains) ) ALLOCATE( bundle%Z(1:bundle%tot_n_domains)) ALLOCATE( bundle%Y(1:bundle%tot_n_domains)) ALLOCATE( bundle%terminal_conductor_list(1:bundle%tot_n_domains)) n_domains=bundle%tot_n_domains do domain=1,bundle%tot_n_domains read(file_unit,*,ERR=9000) ! comment line read(file_unit,*)bundle%n_conductors(domain) read(file_unit,*,ERR=9000) ! comment line read(file_unit,*,ERR=9000)bundle%terminal_conductor_list(domain)%n_elements ALLOCATE( bundle%terminal_conductor_list(domain)%element(1:bundle%terminal_conductor_list(domain)%n_elements) ) do conductor=1,bundle%terminal_conductor_list(domain)%n_elements read(file_unit,*,ERR=9000)bundle%terminal_conductor_list(domain)%element(conductor) end do read(file_unit,*,ERR=9000) ! comment line read(file_unit,*)matrix_dimension bundle%L(domain)%dim=matrix_dimension ALLOCATE( bundle%L(domain)%mat(1:matrix_dimension,1:matrix_dimension) ) CALL dread_matrix(bundle%L(domain)%mat,matrix_dimension,matrix_dimension,matrix_dimension,file_unit) write(*,*)'L=' CALL dwrite_matrix(bundle%L(domain)%mat,matrix_dimension,matrix_dimension,matrix_dimension,0) ! ****** read(file_unit,*,ERR=9000) ! comment line read(file_unit,*)matrix_dimension bundle%C(domain)%dim=matrix_dimension ALLOCATE( bundle%C(domain)%mat(1:matrix_dimension,1:matrix_dimension) ) CALL dread_matrix(bundle%C(domain)%mat,matrix_dimension,matrix_dimension,matrix_dimension,file_unit) write(*,*)'C=' CALL dwrite_matrix(bundle%C(domain)%mat,matrix_dimension,matrix_dimension,matrix_dimension,0) ! ****** read(file_unit,*,ERR=9000) ! comment line CALL read_Sfilter_matrix( bundle%Z(domain),file_unit ) read(file_unit,*,ERR=9000) ! comment line CALL read_Sfilter_matrix( bundle%Y(domain),file_unit ) end do read(file_unit,*,ERR=9000) ! comment line read(file_unit,*)matrix_dimension bundle%global_MI%dim=matrix_dimension ALLOCATE( bundle%global_MI%mat(1:matrix_dimension,1:matrix_dimension) ) CALL dread_matrix(bundle%global_MI%mat,matrix_dimension,matrix_dimension,matrix_dimension,file_unit) read(file_unit,*,ERR=9000) ! comment line read(file_unit,*)matrix_dimension bundle%global_MV%dim=matrix_dimension ALLOCATE( bundle%global_MV%mat(1:matrix_dimension,1:matrix_dimension) ) CALL dread_matrix(bundle%global_MV%mat,matrix_dimension,matrix_dimension,matrix_dimension,file_unit) read(file_unit,*,ERR=9000) ! comment line read(file_unit,*)matrix_dimension bundle%global_L%dim=matrix_dimension ALLOCATE( bundle%global_L%mat(1:matrix_dimension,1:matrix_dimension) ) CALL dread_matrix(bundle%global_L%mat,matrix_dimension,matrix_dimension,matrix_dimension,file_unit) read(file_unit,*,ERR=9000) ! comment line read(file_unit,*)matrix_dimension bundle%global_C%dim=matrix_dimension ALLOCATE( bundle%global_C%mat(1:matrix_dimension,1:matrix_dimension) ) CALL dread_matrix(bundle%global_C%mat,matrix_dimension,matrix_dimension,matrix_dimension,file_unit) read(file_unit,*,ERR=9000) ! comment line CALL read_Sfilter_matrix( bundle%global_Z,file_unit ) read(file_unit,*,ERR=9000) ! comment line CALL read_Sfilter_matrix( bundle%global_Y,file_unit ) ! read the loss model information for each conductor ALLOCATE( bundle%conductor_impedance(1:bundle%tot_n_conductors) ) read(file_unit,*) ! comment line do conductor=1,bundle%tot_n_conductors CALL read_conductor_impedance_model(bundle%conductor_impedance(conductor),file_unit) end do ! read the x and y position for each conductor: this is needed for incident field excitation ALLOCATE( bundle%conductor_x_offset(1:bundle%tot_n_conductors) ) ALLOCATE( bundle%conductor_y_offset(1:bundle%tot_n_conductors) ) read(file_unit,*) ! comment line do conductor=1,bundle%tot_n_conductors read(file_unit,*)bundle%conductor_x_offset(conductor),bundle%conductor_y_offset(conductor) end do ! numbering information read(file_unit,*) ! comment line ALLOCATE( bundle%terminal_conductor_is_shield_flag(1:bundle%tot_n_conductors) ) ALLOCATE( bundle%terminal_conductor_to_inner_domain(1:bundle%tot_n_conductors) ) ALLOCATE( bundle%terminal_conductor_to_outer_domain(1:bundle%tot_n_conductors) ) ALLOCATE( bundle%terminal_conductor_to_global_domain_conductor(1:bundle%tot_n_conductors) ) ALLOCATE( bundle%terminal_conductor_to_local_domain_conductor(1:bundle%tot_n_conductors) ) ALLOCATE( bundle%terminal_conductor_to_reference_terminal_conductor(1:bundle%tot_n_conductors) ) do conductor=1,bundle%tot_n_conductors read(file_unit,*)bundle%terminal_conductor_is_shield_flag(conductor), & bundle%terminal_conductor_to_inner_domain(conductor), & bundle%terminal_conductor_to_outer_domain(conductor), & bundle%terminal_conductor_to_global_domain_conductor(conductor), & bundle%terminal_conductor_to_local_domain_conductor(conductor), & bundle%terminal_conductor_to_reference_terminal_conductor(conductor) end do ! read the conductor labels ALLOCATE( bundle%conductor_label(1:bundle%tot_n_conductors) ) read(file_unit,*) ! comment line do conductor=1,bundle%tot_n_conductors read(file_unit,'(A)')bundle%conductor_label(conductor) end do ! close .bundle file CLOSE(unit=file_unit) if (verbose) write(*,*)'Closed file:',trim(filename) RETURN 9000 run_status='ERROR reading the bundle file:'//trim(filename) CALL write_program_status() STOP 1 END SUBROUTINE read_cable_bundle ! NAME ! SUBROUTINE write_cable_bundle ! ! write the cable bundle structure to a specified unit ! ! COMMENTS ! ! ! HISTORY ! started 2/12/2015 CJS ! stage 3 developments started 21/3/2015 CJS ! generalise conductor impedance for transfer impedance model (stage 5) CJS 12/5/2016 ! SUBROUTINE write_cable_bundle(bundle,file_unit) USE type_specifications USE constants USE general_module USE cable_module USE maths IMPLICIT NONE ! variables passed to subroutine type(bundle_specification_type),intent(INOUT) :: bundle integer,intent(IN) :: file_unit ! local variables character(len=filename_length) :: filename character(len=line_length) :: line integer :: cable integer :: domain integer :: matrix_dimension integer :: conductor integer :: n_cables_to_write ! START ! Open the output (.bundle) file filename=trim(MOD_bundle_lib_dir)//trim(bundle%bundle_name)//bundle_file_extn open(unit=file_unit,file=trim(filename)) ! write .bundle file write(file_unit,'(A)')trim(bundle%version) write(file_unit,'(A)')trim(bundle%bundle_name) write(file_unit,*)bundle%n_cables_without_ground_plane,' ! number of cables not including ground plane' write(file_unit,*)bundle%n_cables,' ! number of cables, cable name and x y coordinates follow...' ! write cable names and coordinates do cable=1,bundle%n_cables_without_ground_plane write(file_unit,'(A)')trim(bundle%cable(cable)%cable_name) write(file_unit,*)bundle%cable_x_offset(cable),bundle%cable_y_offset(cable), & bundle%cable_angle(cable),' x y coordinates and angle of cable ' end do ! next cable ! write ground plane specification if (bundle%ground_plane_present) then write(file_unit,'(A)')'ground_plane' ! note: convert angle back to degrees write(file_unit,*)bundle%ground_plane_x,bundle%ground_plane_y,bundle%ground_plane_angle*180d0/pi, & ' x y coordinates and angle of ground plane ' write(file_unit,*)bundle%ground_plane_nx,bundle%ground_plane_ny,' ground plane normal direction' write(file_unit,*)bundle%ground_plane_cable_side,' orientation of cables wrt ground plane' write(file_unit,*)bundle%ground_plane_sigma,' ground plane conductivity' write(file_unit,*)bundle%ground_plane_w,' ground plane width ' write(file_unit,*)bundle%ground_plane_h,' ground plane height (thickness) ' write(file_unit,*)bundle%ground_plane_Rdc,' ground plane Rdc ' write(file_unit,*)bundle%ground_plane_nsegx,' ground plane nsegx ' write(file_unit,*)bundle%ground_plane_nsegz,' ground plane nsegz ' write(file_unit,*)bundle%ground_plane_nh,' ground plane nh ' else write(file_unit,'(A)')'no_ground_plane' end if write(file_unit,*)bundle%tot_n_conductors,' # total number of conductors' write(file_unit,*)bundle%tot_n_external_conductors,' # total number of external conductors' write(file_unit,*)bundle%system_dimension,' # dimension of the matrix system characterising the MTL propagation' write(file_unit,*)bundle%tot_n_domains,' # number of domains' do domain=1,bundle%tot_n_domains write(file_unit,*)'Domain number',domain write(file_unit,*)bundle%n_conductors(domain),' number of conductors in this domain' write(file_unit,*)'terminal_conductor_list' write(file_unit,*)bundle%terminal_conductor_list(domain)%n_elements,' ! number of elements ' do conductor=1,bundle%terminal_conductor_list(domain)%n_elements write(file_unit,*)bundle%terminal_conductor_list(domain)%element(conductor) end do matrix_dimension=bundle%L(domain)%dim write(file_unit,*)'Per-Unit-length Inductance Matrix, [L]' write(file_unit,*)matrix_dimension,' Dimension of [L]' CALL dwrite_matrix(bundle%L(domain)%mat,matrix_dimension,matrix_dimension,matrix_dimension,file_unit) write(file_unit,*)'Per-Unit-length Capacitance Matrix, [C]' write(file_unit,*)matrix_dimension,' Dimension of [C]' CALL dwrite_matrix(bundle%C(domain)%mat,matrix_dimension,matrix_dimension,matrix_dimension,file_unit) write(file_unit,*)'Per-Unit-length Impedance Matrix, Z' CALL write_Sfilter_matrix( bundle%Z(domain),file_unit ) write(file_unit,*)'Per-Unit-length Admittance Matrix, Y' CALL write_Sfilter_matrix( bundle%Y(domain),file_unit ) end do write(file_unit,*)'# Global current domain transformation matrix, [MI]' matrix_dimension=bundle%global_MI%dim-1 write(file_unit,*)matrix_dimension,' Dimension of [MI]' CALL dwrite_matrix(bundle%global_MI%mat,matrix_dimension,matrix_dimension,bundle%global_MI%dim,file_unit) write(file_unit,*)'# Global voltage domain transformation matrix, [MV]' matrix_dimension=bundle%global_MV%dim-1 write(file_unit,*)matrix_dimension,' Dimension of [MV]' CALL dwrite_matrix(bundle%global_MV%mat,matrix_dimension,matrix_dimension,bundle%global_MI%dim,file_unit) write(file_unit,*)'# Global domain based inductance matrix, [L]' matrix_dimension=bundle%global_L%dim write(file_unit,*)matrix_dimension,' Dimension of [L]' CALL dwrite_matrix(bundle%global_L%mat,matrix_dimension,matrix_dimension,matrix_dimension,file_unit) write(file_unit,*)'# Global domain based capacitance matrix, [C]' matrix_dimension=bundle%global_C%dim write(file_unit,*)matrix_dimension,' Dimension of [C]' CALL dwrite_matrix(bundle%global_C%mat,matrix_dimension,matrix_dimension,matrix_dimension,file_unit) write(file_unit,*)'Global domain based Per-Unit-length Impedance Matrix, Z' CALL write_Sfilter_matrix( bundle%global_Z,file_unit ) write(file_unit,*)'Global domain based Per-Unit-length Admittance Matrix, Y' CALL write_Sfilter_matrix( bundle%global_Y,file_unit ) ! write the loss model information for each conductor write(file_unit,*)' # conductor impedance models' do conductor=1,bundle%tot_n_conductors CALL write_conductor_impedance_model(bundle%conductor_impedance(conductor),file_unit) end do ! write the x and y position for each conductor: this is needed for incident field excitation write(file_unit,*)' # conductor x y positions' do conductor=1,bundle%tot_n_conductors write(file_unit,*)bundle%conductor_x_offset(conductor),bundle%conductor_y_offset(conductor) end do ! numbering information write(file_unit,'(A)')'is_shield tc_to_in_domain tc_to_out_domain tc_to_gdc tc_to_ldc tc_to_ref_tc ' do conductor=1,bundle%tot_n_conductors write(file_unit,8000)bundle%terminal_conductor_is_shield_flag(conductor), & bundle%terminal_conductor_to_inner_domain(conductor), & bundle%terminal_conductor_to_outer_domain(conductor), & bundle%terminal_conductor_to_global_domain_conductor(conductor), & bundle%terminal_conductor_to_local_domain_conductor(conductor), & bundle%terminal_conductor_to_reference_terminal_conductor(conductor) 8000 format(L5,5I15) end do ! write the conductor labels write(file_unit,*)' # Conductor labels' do conductor=1,bundle%tot_n_conductors write(file_unit,*)trim(bundle%conductor_label(conductor)) end do ! close .bundle file CLOSE(unit=file_unit) RETURN END SUBROUTINE write_cable_bundle ! NAME ! SUBROUTINE deallocate_cable_bundle ! ! deallocate cable_bundle structure data ! ! COMMENTS ! ! ! HISTORY ! started 2/12/2015 CJS ! generalise conductor impedance for transfer impedance model (stage 5) CJS 12/5/2016 ! SUBROUTINE deallocate_cable_bundle(bundle) USE type_specifications USE general_module USE cable_module IMPLICIT NONE ! variables passed to subroutine type(bundle_specification_type),intent(INOUT) :: bundle ! local variables integer :: cable integer :: domain integer :: conductor ! START do cable=1,bundle%n_cables CALL deallocate_cable(bundle%cable(cable)) end do if (allocated(bundle%cable)) DEALLOCATE( bundle%cable ) if (allocated(bundle%cable_x_offset )) DEALLOCATE( bundle%cable_x_offset ) if (allocated(bundle%cable_y_offset )) DEALLOCATE( bundle%cable_y_offset ) if (allocated(bundle%cable_angle )) DEALLOCATE( bundle%cable_angle ) if (allocated(bundle%n_conductors )) DEALLOCATE( bundle%n_conductors ) if (allocated(bundle%terminal_conductor_list )) then do domain=1,bundle%tot_n_domains if (allocated(bundle%terminal_conductor_list(domain)%element )) & DEALLOCATE( bundle%terminal_conductor_list(domain)%element ) end do DEALLOCATE( bundle%terminal_conductor_list ) end if if (allocated(bundle%L )) then do domain=1,bundle%tot_n_domains if (allocated(bundle%L(domain)%mat )) DEALLOCATE( bundle%L(domain)%mat ) end do DEALLOCATE( bundle%L ) end if if (allocated(bundle%C )) then do domain=1,bundle%tot_n_domains if (allocated(bundle%C(domain)%mat )) DEALLOCATE( bundle%C(domain)%mat ) end do DEALLOCATE( bundle%C ) end if if (ALLOCATED(bundle%Z)) then do domain=1,bundle%tot_n_domains CALL deallocate_Sfilter_matrix( bundle%Z(domain) ) end do DEALLOCATE(bundle%Z) end if if (ALLOCATED(bundle%Y)) then do domain=1,bundle%tot_n_domains CALL deallocate_Sfilter_matrix( bundle%Y(domain) ) end do DEALLOCATE(bundle%Y) end if if (allocated(bundle%global_MI%mat)) DEALLOCATE( bundle%global_MI%mat ) if (allocated(bundle%global_MV%mat)) DEALLOCATE( bundle%global_MV%mat ) if (allocated(bundle%global_L%mat)) DEALLOCATE( bundle%global_L%mat ) if (allocated(bundle%global_C%mat)) DEALLOCATE( bundle%global_C%mat ) CALL deallocate_Sfilter_matrix( bundle%global_Z ) CALL deallocate_Sfilter_matrix( bundle%global_Y ) if (ALLOCATED(bundle%conductor_impedance)) then do conductor=1,bundle%tot_n_conductors CALL deallocate_conductor_impedance_model(bundle%conductor_impedance(conductor)) end do DEALLOCATE(bundle%conductor_impedance) end if if (ALLOCATED( bundle%conductor_x_offset )) DEALLOCATE( bundle%conductor_x_offset ) if (ALLOCATED( bundle%conductor_y_offset )) DEALLOCATE( bundle%conductor_y_offset ) ! Numbering information required for transfer impedance calculation if(ALLOCATED( bundle%terminal_conductor_is_shield_flag )) DEALLOCATE( bundle%terminal_conductor_is_shield_flag ) if(ALLOCATED( bundle%terminal_conductor_to_inner_domain )) DEALLOCATE( bundle%terminal_conductor_to_inner_domain ) if(ALLOCATED( bundle%terminal_conductor_to_outer_domain )) DEALLOCATE( bundle%terminal_conductor_to_outer_domain ) if(ALLOCATED( bundle%terminal_conductor_to_global_domain_conductor )) & DEALLOCATE( bundle%terminal_conductor_to_global_domain_conductor ) if(ALLOCATED( bundle%terminal_conductor_to_local_domain_conductor )) & DEALLOCATE( bundle%terminal_conductor_to_local_domain_conductor ) if(ALLOCATED( bundle%terminal_conductor_to_reference_terminal_conductor )) & DEALLOCATE( bundle%terminal_conductor_to_reference_terminal_conductor ) if (ALLOCATED( bundle%conductor_label )) then DEALLOCATE(bundle%conductor_label) end if RETURN END SUBROUTINE deallocate_cable_bundle ! ! NAME ! SUBROUTINE check_cables_wrt_ground_plane ! ! check that the cables are all on one side of the ground plane to ensure consistency ! and indicate which side the cables are on by setting bundle%ground_plane_cable_side ! ! COMMENTS ! ! ! HISTORY ! started 29/06/2016 CJS ! SUBROUTINE check_cables_wrt_ground_plane(bundle) USE type_specifications USE general_module USE cable_module IMPLICIT NONE ! variables passed to subroutine type(bundle_specification_type),intent(INOUT) :: bundle ! local variables integer :: cable real(dp) :: norm_x,norm_y real(dp) :: p real(dp) :: gp_offset real(dp) :: gp_offset_cable_1 ! START ! calculate the ground plane normal direction norm_x=bundle%ground_plane_nx norm_y=bundle%ground_plane_ny ! The equation of the ground plane is r.norm-P=0, calculate P here P=bundle%ground_plane_x*norm_x+bundle%ground_plane_y*norm_y ! loop over all of the cables, calculating the offset from the ground plane in the normal direction ! offset=r.norm-P, and check that all offsets have the same sign if (verbose) then write(*,*)'Checking the orientation of cables wrt the ground plane' write(*,*)'GP normal:',norm_x,norm_y write(*,*)'GP offset:',bundle%ground_plane_x+bundle%ground_plane_y write(*,*)'P=',P end if do cable=1,bundle%n_cables_without_ground_plane gp_offset=bundle%cable_x_offset(cable)*norm_x+bundle%cable_y_offset(cable)*norm_y-P if (verbose) then write(*,*)'Cable',cable,' offset',gp_offset end if if (cable.EQ.1) then gp_offset_cable_1=gp_offset else ! check consistency of ground plane offsets if (gp_offset_cable_1*gp_offset.LE.0d0) then run_status='ERROR there are cables both sides of the ground plane' CALL write_program_status() STOP 1 end if end if ! cable.NE.1 end do ! next cable to check if (gp_offset_cable_1.GT.0d0) then bundle%ground_plane_cable_side=1 else bundle%ground_plane_cable_side=-1 end if RETURN END SUBROUTINE check_cables_wrt_ground_plane ! ! NAME ! SUBROUTINE check_cable_intersections ! ! check that the cables do not intersect ! the intersection test checks intersection of the outer shape of each cable ! which is usually the outer dielectric layer apart from overshields or D connectors which ! are assumed not to have dielectric ! ! COMMENTS ! Revised to work with all cable types ! ! HISTORY ! started 8/10/2016 CJS ! ! SUBROUTINE check_cable_intersection(bundle) USE type_specifications USE general_module USE cable_module IMPLICIT NONE ! variables passed to subroutine type(bundle_specification_type),intent(IN) :: bundle ! local variables integer :: cable1,cable2 integer :: nec1,nec2,ec1,ec2 integer :: shape1,shape2 integer :: type1,type2 real(dp) :: ox,oy,theta real(dp) :: r real(dp) :: wd,hd integer :: nc,ncrow(2) real(dp) :: rw,p,s,o,W(2) integer :: npts1 real(dp),allocatable :: shape1_x(:) real(dp),allocatable :: shape1_y(:) integer :: npts2 real(dp),allocatable :: shape2_x(:) real(dp),allocatable :: shape2_y(:) logical :: intersect logical :: intersect2 logical :: nested_1_in_2 logical :: nested_2_in_1 logical :: gp_intersect logical :: intersection_found ! START intersection_found=.FALSE. ! loop over cable 1 do cable1=1,bundle%n_cables_without_ground_plane ox=bundle%cable_x_offset(cable1) oy=bundle%cable_y_offset(cable1) theta=bundle%cable_angle(cable1) type1=bundle%cable(cable1)%cable_type if ((bundle%cable(cable1)%cable_type.NE.cable_geometry_type_flex_cable).AND. & (bundle%cable(cable1)%cable_type.NE.cable_geometry_type_ML_flex_cable) ) then nec1=bundle%cable(cable1)%n_external_conductors else nec1=1 end if ! loop over the external conductors of cable 1 do ec1=1,nec1 shape1=bundle%cable(cable1)%external_model(ec1)%conductor_type ! generate a list of points on cable 1 outer surface if (shape1.EQ.rectangle) then wd=bundle%cable(cable1)%external_model(1)%dielectric_width hd=bundle%cable(cable1)%external_model(1)%dielectric_height CALL generate_rectangle_points(npts1,shape1_x,shape1_y,ox,oy,theta,wd,hd) else if (shape1.EQ.circle) then r=bundle%cable(cable1)%external_model(ec1)%dielectric_radius CALL generate_circle_points(npts1,shape1_x,shape1_y,ox,oy,r) else if (shape1.EQ.Dshape) then nc=bundle%cable(cable1)%tot_n_conductors rw=bundle%cable(cable1)%parameters(1) p=bundle%cable(cable1)%parameters(2) s=bundle%cable(cable1)%parameters(3) o=bundle%cable(cable1)%parameters(4) ncrow(1)=nc/2 ncrow(2)=(nc-1)-ncrow(1) w(1)=(ncrow(1)-1)*p w(2)=(ncrow(2)-1)*p CALL generate_Dshape_points(npts1,shape1_x,shape1_y,ox,oy,w(1)/2d0,w(2)/2d0,s/2d0,rw+o,theta) end if if (bundle%ground_plane_present) then ! check for intersection of the cables with the ground plane gp_intersect=.FALSE. if (verbose) write(*,*)'TESTING: intersection of cable',cable1,' with the ground plane' CALL gptest(npts1,shape1_x,shape1_y,gp_intersect) if (gp_intersect) then if (verbose) write(*,*)'Cable ',cable1,' intersects the ground plane' intersection_found=.TRUE. end if end if ! loop over cable 2 do cable2=cable1+1,bundle%n_cables_without_ground_plane ox=bundle%cable_x_offset(cable2) oy=bundle%cable_y_offset(cable2) theta=bundle%cable_angle(cable2) type2=bundle%cable(cable2)%cable_type if ((bundle%cable(cable1)%cable_type.NE.cable_geometry_type_flex_cable).AND. & (bundle%cable(cable1)%cable_type.NE.cable_geometry_type_ML_flex_cable) ) then nec2=bundle%cable(cable2)%n_external_conductors else nec2=1 end if ! loop over the external conductors of cable 2 do ec2=1,nec2 shape2=bundle%cable(cable2)%external_model(ec2)%conductor_type ! generate a list of points on cable 1 outer surface if (shape2.EQ.rectangle) then wd=bundle%cable(cable2)%external_model(1)%dielectric_width hd=bundle%cable(cable2)%external_model(1)%dielectric_height CALL generate_rectangle_points(npts2,shape2_x,shape2_y,ox,oy,theta,wd,hd) else if (shape2.EQ.circle) then r=bundle%cable(cable2)%external_model(ec2)%dielectric_radius CALL generate_circle_points(npts2,shape2_x,shape2_y,ox,oy,r) else if (shape2.EQ.Dshape) then nc=bundle%cable(cable2)%tot_n_conductors rw=bundle%cable(cable2)%parameters(1) p=bundle%cable(cable2)%parameters(2) s=bundle%cable(cable2)%parameters(3) o=bundle%cable(cable2)%parameters(4) ncrow(1)=nc/2 ncrow(2)=(nc-1)-ncrow(1) w(1)=(ncrow(1)-1)*p w(2)=(ncrow(2)-1)*p CALL generate_Dshape_points(npts2,shape2_x,shape2_y,ox,oy,w(1)/2d0,w(2)/2d0,s/2d0,rw+o,theta) end if ! assume that there is no intersection and the shapes are not nested intersect=.FALSE. nested_1_in_2=.FALSE. nested_2_in_1=.FALSE. ! check for the intersection of cable 1 with cable 2 and flag if cable 1 is nested within cable 2. if (verbose) write(*,*)'TESTING: intersection of cable',cable1,' with cable ',cable2 CALL shape_test(npts1,shape1_x,shape1_y,npts2,shape2_x,shape2_y,intersect,nested_1_in_2) if (intersect) then if (verbose) write(*,*)'Cable ',cable1,' intersects Cable',cable2 intersection_found=.TRUE. end if ! Cables can only be nested if the outer cable is an overshield if ( nested_1_in_2 ) then if (type1.NE.cable_geometry_type_overshield) then if (verbose) write(*,*)'Cable ',cable2,' is inside Cable',cable1 intersection_found=.TRUE. else if (verbose) write(*,*)'Cable ',cable2,' is OK inside Overshield Cable',cable1 end if end if ! check for the intersection of cable 1 with cable 2 and flag if cable 1 is nested within cable 2. if (verbose) write(*,*)'TESTING: intersection of cable',cable2,' with cable ',cable1 CALL shape_test(npts2,shape2_x,shape2_y,npts1,shape1_x,shape1_y,intersect2,nested_2_in_1) if (intersect) then if (verbose) write(*,*)'Cable ',cable2,' intersects Cable',cable1 intersection_found=.TRUE. end if ! Cables can only be nested if the outer cable is an overshield if ( nested_2_in_1 ) then if (type2.NE.cable_geometry_type_overshield) then if (verbose) write(*,*)'Cable ',cable1,' is inside Cable',cable2 intersection_found=.TRUE. else if (verbose) write(*,*)'Cable ',cable1,' is OK inside Overshield Cable',cable2 end if end if DEALLOCATE( shape2_x ) DEALLOCATE( shape2_y ) end do ! next external conductor of cable 2 end do ! next cable 2 DEALLOCATE( shape1_x ) DEALLOCATE( shape1_y ) end do ! next external conductor of cable1 end do ! next cable1 if (intersection_found) then run_status='ERROR there are cables in the bundle which intersect' CALL write_program_status() STOP 1 end if RETURN END SUBROUTINE check_cable_intersection ! ! NAME ! SUBROUTINE gptest ! ! check that a cable does not intersect the ground plane ! ! COMMENTS ! ! ! HISTORY ! started 20/4/2017 CJS ! ! SUBROUTINE gptest(npts,shape_x,shape_y,gp_intersect) USE type_specifications IMPLICIT NONE integer,intent(IN) :: npts real(dp),intent(IN) :: shape_x(npts) real(dp),intent(IN) :: shape_y(npts) logical,intent(OUT) :: gp_intersect ! local variables integer :: i ! START do i=1,npts if (shape_y(i).LE.0D0) then gp_intersect=.TRUE. RETURN end if end do ! next point RETURN END SUBROUTINE gptest ! ! NAME ! SUBROUTINE shape_test ! ! check that cables do not intersect, also flag whether the first cable is nested ! i.e. lies entirely within the second ! ! COMMENTS ! We assume that the shape of the cable is convex and that the shape is closed ! ! HISTORY ! started 20/4/2017 CJS ! ! SUBROUTINE shape_test(npts1,shape1_x,shape1_y,npts2,shape2_x,shape2_y,intersect,nested) USE type_specifications IMPLICIT NONE integer,intent(IN) :: npts1 real(dp),intent(IN) :: shape1_x(npts1) real(dp),intent(IN) :: shape1_y(npts1) integer,intent(IN) :: npts2 real(dp),intent(IN) :: shape2_x(npts2) real(dp),intent(IN) :: shape2_y(npts2) logical,intent(OUT) :: intersect ! flag to indicate that shapes 1 and 2 intersect logical,intent(OUT) :: nested ! flag to indicate shape 2 lies within shape 1 ! local variables integer :: i,j real(dp) :: vx1,vy1 real(dp) :: vx2,vy2 real(dp) :: vp real(dp) :: dirn logical :: inside_point_found logical :: outside_point_found logical :: local_outside_point_found ! START intersect=.FALSE. nested=.FALSE. inside_point_found=.FALSE. outside_point_found=.FALSE. ! loop around test points in shape 2 do j=1,npts2 local_outside_point_found=.FALSE. dirn=0.0 ! set to zero initially to indicate that it is not known ! loop around line segments in shape 1. Note there is one less line segment than the number of points do i=1,npts1-1 ! vector along line segment in shape 1 vx1=shape1_x(i+1)-shape1_x(i) vy1=shape1_y(i+1)-shape1_y(i) ! vector from first point in line segment to test point vx2=shape2_x(j)-shape1_x(i) vy2=shape2_y(j)-shape1_y(i) ! calculate the z component of the vector product of the two vectors vp=vx1*vy2-vx2*vy1 if (dirn.EQ.0D0) then ! maybe we need to use abs(dirn).LT.small here... ! this is the first value so this sets the direction to test for all other sets of points dirn=vp else ! test whether the sign of the direction has changed if (dirn*vp.LT.0.0) then ! the signs are different then this indicates a test point outside shape 1 local_outside_point_found=.TRUE. end if ! change in sign of direction end if ! first test line segment. end do ! next line segment in shape 1 if (local_outside_point_found) then outside_point_found=.TRUE. else inside_point_found=.TRUE. end if end do ! next test point in shape 2 ! work out whether we have an intersection if (inside_point_found.AND.outside_point_found) intersect=.TRUE. ! work out whether we have a nested shape if (inside_point_found.AND.(.NOT.outside_point_found)) nested=.TRUE. RETURN END SUBROUTINE END MODULE cable_bundle_module