module envelope_beambeam_setup_mod contains !+ ! subroutine ENVELOPE_BEAMBEAM_SETUP (ring, particle, current, bunchwt, coupling, slices) ! ! Subroutine to add beambeam element at IP ! ! Input: ! RING -- Lat_struct: Ring containing the lattice to be modified ! particle -- Integer : Strong beam particle type (+1:e+, -1:e-) ! current -- real : bunch current in ma ! bunchwt -- real : Weight for the bunch current of the colliding bunch ! This allows it to be different from the parasitic bbi ! A value of unity gives the full CURRENT value. ! Since BUNCHWT is used to multiply VALUE(CHARGE$) ! of the IP BBI element, which is the only ! element-specific parameter used to weight the strength ! of the interaction, it permits setting a strength ! which differs from the strengths of the parasitic BBI. ! 25 Aug 2004 jac ! slices -- integer, optional : number of slices, default is 1 ! coupling -- real: xy coupling of strong beam emittance ! Output: ! RING -- Lat_struct: Ring includes Beam-beam element at IP ! subroutine envelope_beambeam_setup(ring, particle, current, bunchwt, coupling, slices) use bmad use bmadz_struct use bmadz_interface use superimpose_mod implicit none type ( lat_struct ) ring type (lat_struct), save :: ring_oppos type (coord_struct), allocatable :: co_(:), co_oppos_(:) type (ele_struct) ele, beambeam_ele type (normal_modes_struct) mode integer particle integer ierr integer i,n,ix_ip integer, optional, intent(in) :: slices logical err_flag real(rp) current, bunchwt real(rp) major, minor, theta, vert_size_coupled, vert_size_emit real(rp) coupling allocate( co_(0:ring%n_ele_max)) allocate( co_oppos_(0:ring%n_ele_max)) ! init ring%param%particle = particle call twiss_at_start(ring) ! type *,' envelope_beambeam_setup:1 beta ',ring%ele(0)%a%beta, ring%ele(0)%b%beta co_(0)%vec = 0. ! call closed_orbit_at_start(ring, co_(0), 4, .true.) ! type *, ' closed_orbit: co_',co_(0)%vec(1:4) ! call track_all (ring, co_) call closed_orbit_calc(ring, co_, 4) call lat_make_mat6(ring,-1,co_) call twiss_at_start(ring) ! type *,' envelope_beambeam_setup:2 beta ',ring%ele(0)%a%beta, ring%ele(0)%b%beta call twiss_propagate_all (ring) ! set up for lrbbi if( current == 0)then type *,' ENVELOPE_BEAMBEAM_SETUP: bunch specifications incomplete ' stop endif ring%param%n_part = current*0.001 *(ring%param%total_length/c_light)/e_charge ring_oppos = ring ring_oppos%param%n_part = 0. ring_oppos%param%particle = -particle call lat_make_mat6(ring_oppos,-1) co_oppos_(0)%vec = 0. ! call closed_orbit_at_start(ring_oppos, co_oppos_(0), 4, .true.) call closed_orbit_calc(ring_oppos, co_oppos_, 4) call track_all (ring_oppos, co_oppos_) call lat_make_mat6(ring_oppos,-1,co_oppos_) call twiss_at_start(ring_oppos) ! type *,' envelope_beambeam_setup:3 beta ',ring%ele(0)%a%beta, ring%ele(0)%b%beta call twiss_propagate_all (ring_oppos) call set_on_off (rfcavity$, ring_oppos, on$) call set_z_tune(ring_oppos%branch(0), ring_oppos%z%tune) call radiation_integrals (ring_oppos, co_oppos_, mode) call ellipse (ring_oppos%ele(0), major, minor, theta) call init_ele (beambeam_ele) beambeam_ele%name = 'IP_COLLISION' beambeam_ele%key = beambeam$ ! beambeam_ele%value(sig_x$) = sqrt(ring_oppos%ele(0)%a%beta * mode%a%emittance & ! + (ring_oppos%ele(0)%a%eta * mode%sige_e)**2) beambeam_ele%value(sig_x$) = sqrt( mode%a%emittance * major**2 + & (ring_oppos%ele(0)%a%eta * mode%sige_e)**2) ! beambeam_ele%value(sig_y$) = sqrt(ring_oppos%ele(0)%b%beta * mode%a%emittance*0.01 & ! + (ring_oppos%ele(0)%b%eta * mode%sige_e)**2) vert_size_coupled = mode%a%emittance * max( minor**2, coupling * ring_oppos%ele(0)%b%beta) vert_size_emit = mode%b%emittance * ring_oppos%ele(0)%b%beta beambeam_ele%value(sig_y$) = sqrt(vert_size_coupled + vert_size_emit & + (ring_oppos%ele(0)%b%eta * mode%sige_e)**2) beambeam_ele%value(sig_z$) = mode%sig_z if(present(slices))then beambeam_ele%value(n_slice$) = slices else beambeam_ele%value(n_slice$) = 1 endif beambeam_ele%value(charge$) = -1 * bunchwt beambeam_ele%value(x_pitch$) = co_oppos_(0)%vec(2) beambeam_ele%value(y_pitch$) = co_oppos_(0)%vec(4) beambeam_ele%value(x_offset$) = co_oppos_(0)%vec(1) beambeam_ele%value(y_offset$) = co_oppos_(0)%vec(3) beambeam_ele%value(tilt$) = theta !tilt of beam ellipse due to coupling ix_ip = 1 call add_superimpose (ring, beambeam_ele, ix_ip, err_flag) call lat_make_mat6(ring, -1) type '(1x,a14)', ' Strong beam: ' type '(1x,a12,e12.4)', ' sigma_x = ',beambeam_ele%value(sig_x$) type '(1x,a12,e12.4)', ' sigma_y = ',beambeam_ele%value(sig_y$) type '(1x,a12,e12.4)', ' sigma_z = ',beambeam_ele%value(sig_z$) type '(1x,a14,e12.4,a4,e12.4)', ' Pitch : x= ',beambeam_ele%value(x_pitch$), & ' y= ',beambeam_ele%value(y_pitch$) type '(1x,a14,e12.4,a4,e12.4)', ' Offset : x= ',beambeam_ele%value(x_offset$), & ' y= ',beambeam_ele%value(y_offset$) type '(1x,a9,e12.4)', ' Tilt = ', beambeam_ele%value(tilt$) type '(1x,a11,e12.4)', ' Charge = ', beambeam_ele%value(charge$) ! type '(1x,a1,5x,2x,a12,7x,a1,4x,4a9)','n','Element name','s',' x_pitch ',' y_pitch ', & ! 'x_offset ', & ! 'y_offset' ! do i=1, ring%n_ele_track ! ele = ring%ele(i) ! if(ele%key == beambeam$)then ! n=n+1 ! type '(i3,4x,a16,4x,f7.3,2x,4f9.5)',n,ring%ele(i)%name, ring%ele(i)%s, & ! ele%value(x_pitch$), ele%value(y_pitch$), ele%value(x_offset$), ele%value(y_offset$) ! endif ! end do deallocate( co_) deallocate( co_oppos_) end subroutine end module