diff --git a/bmad/code/pointer_to_attribute.f90 b/bmad/code/pointer_to_attribute.f90 index 7edc09d8d6..9da6fd0cfe 100644 --- a/bmad/code/pointer_to_attribute.f90 +++ b/bmad/code/pointer_to_attribute.f90 @@ -533,26 +533,21 @@ subroutine pointer_to_attribute (ele, attrib_name, do_allocation, a_ptr, err_fla case ('N_LORD'); a_ptr%i => ele%n_lord case ('LR_FREQ_SPREAD', 'LR_SELF_WAKE_ON', 'LR_WAKE%AMP_SCALE', 'LR_WAKE%TIME_SCALE', & 'LR_WAKE%FREQ_SPREAD', 'LR_WAKE%SELF_WAKE_ON', & - 'SR_WAKE%SCALE_WITH_LENGTH', 'SR_WAKE%AMP_SCALE', 'SR_WAKE%Z_SCALE') + 'SR_WAKE%SCALE_WITH_LENGTH', 'SR_WAKE%AMP_SCALE', 'SR_WAKE%Z_SCALE', & + 'SR_WAKE%Z_LONG%SMOOTHING_SIGMA') if (.not. associated(ele%wake)) then if (.not. do_allocation) goto 9100 call init_wake (ele%wake, 0, 0, 0, 0, .true.) endif select case (a_name) - case ('LR_SELF_WAKE_ON', 'LR_WAKE%SELF_WAKE_ON') - a_ptr%l => ele%wake%lr%self_wake_on - case ('LR_WAKE%AMP_SCALE') - a_ptr%r => ele%wake%lr%amp_scale - case ('LR_WAKE%TIME_SCALE') - a_ptr%r => ele%wake%lr%time_scale - case ('LR_FREQ_SPREAD', 'LR_WAKE%FREQ_SPREAD') - a_ptr%r => ele%wake%lr%freq_spread - case ('SR_WAKE%AMP_SCALE') - a_ptr%r => ele%wake%sr%amp_scale - case ('SR_WAKE%Z_SCALE') - a_ptr%r => ele%wake%sr%z_scale - case ('SR_WAKE%SCALE_WITH_LENGTH') - a_ptr%l => ele%wake%sr%scale_with_length + case ('SR_WAKE%Z_LONG%SMOOTHING_SIGMA'); a_ptr%r => ele%wake%sr%z_long%smoothing_sigma + case ('SR_WAKE%AMP_SCALE'); a_ptr%r => ele%wake%sr%amp_scale + case ('SR_WAKE%Z_SCALE'); a_ptr%r => ele%wake%sr%z_scale + case ('SR_WAKE%SCALE_WITH_LENGTH'); a_ptr%l => ele%wake%sr%scale_with_length + case ('LR_SELF_WAKE_ON', 'LR_WAKE%SELF_WAKE_ON'); a_ptr%l => ele%wake%lr%self_wake_on + case ('LR_WAKE%AMP_SCALE'); a_ptr%r => ele%wake%lr%amp_scale + case ('LR_WAKE%TIME_SCALE'); a_ptr%r => ele%wake%lr%time_scale + case ('LR_FREQ_SPREAD', 'LR_WAKE%FREQ_SPREAD'); a_ptr%r => ele%wake%lr%freq_spread end select case ('H_MISALIGN%ACTIVE'); a_ptr%l => ele%photon%h_misalign%active diff --git a/bmad/doc/attributes.tex b/bmad/doc/attributes.tex index 48e6aa71e0..94243b17d3 100644 --- a/bmad/doc/attributes.tex +++ b/bmad/doc/attributes.tex @@ -2930,19 +2930,23 @@ \subsection{Short-Range Wakes} for short-range wakefields are given in \Sref{s:sr.wake.eq}. The \vn{sr_wake} attribute is used to set wakefield parameters. The general form of this attribute is: \begin{example} - sr_wake = \{z_max = , z_scale = , amp_scale = , + sr_wake = \{ + z_max = , z_scale = , amp_scale = , scale_with_length = , longitudinal = \{, , , , \}, ... longitudinal = \{...\}, transverse = \{, , , , , \}, ... - transverse = \{...\} \} + transverse = \{...\} + z_long = \{...\} + \} \end{example} The \vn{sr_wake} structure has optional components \vn{z_max}, \vn{z_scale}, \vn{amp_scale}, and -\vn{scale_with_length} along with zero or more \vn{longitudinal} sub-structures each one specifying -a single longitudinal mode, and zero or more \vn{transverse} sub-structures each one specifying a -single transverse mode. Example: +\vn{scale_with_length}. To specify wakes, zero or more \vn{longitudinal} sub-structures each one specifying +a single longitudinal mode, zero or more \vn{transverse} sub-structures each one specifying a +single transverse mode and zero or one \vn{z_long} sub-structures that can be used to specify the longitudinal +wake as a function of z-position (\sref{s:sr.wake.z.long}). Example: \begin{example} cav9: lcavity, ..., sr_wake = \{z_max = 1.3e-3, scale_with_length = F, longitudinal = \{3.23e14, 1.23e3, 3.62e3, 0.123, none\}, @@ -2974,8 +2978,8 @@ \subsection{Short-Range Wakes} Monopole modes are modes that are independent of transverse position and dipole modes are modes that are linear in the transverse position. -For the \vn{longitudinal} sub-structures, there is a 5\Th -component which gives the transverse position dependence of the wake. Possible values are: +For the \vn{longitudinal} sub-structures, there is a 5\Th component, called +\vn{position_dependence}, gives the transverse position dependence of the wake. Possible values are: \begin{example} none ! No position dependence x_leading ! Linear in the leading particle x-position @@ -3026,6 +3030,49 @@ \subsection{Short-Range Wakes} zero and the transverse wake has no terms independent of the transverse offsets nor terms that depend upon the trailing particle offset. +%----------------------------------------------------------------- +\subsection{Longitudinal Short-Range Wake With a Wake table} +\label{s:sr.wake.z.long} + +The \vn{z_long} substructure of \vn{sr_wake} (\sref{s:sr.wake.syntax}) is an alternative to the +\vn{longitudinal} substructure. With the \vn{z_long} substructure, the wake is specified using a +table of equally spaced points in $z$-position or as a function of time. The syntax for \vn{z_long} +is +\begin{example} + z_long = \{ + time_based = , + smoothing_sigma = , + position_dependence = , + w = \{ + , ! Wake table as a function of z-position (or time). + , + ... + \} + \} +\end{example} +The \vn{w} component gives the single particle wake as a function of either z-position or time +depending upon the setting of the optional \vn{time_based} logical. The default is \vn{False}. The z +(or time) values in the \vn{w} table must be equally spaced and in increasing order. Unlike the +\vn{longitudinal} and \vn{transverse} pseudo-mode description, the wake can be finite for both +positive and negative z (or time). When the wake is parsed, if \vn{time_based} is \vn{True}, the +table is converted to be z-based using the conversion $z = -c \, t$ where $c$ is the speed of +light. If the table is not symmetric around $z = 0$ (if there are more points on one size of zero +than the other), the table is extended to be symmetric (this is needed for tracking). During +tracking, the half-width of the wake (the magnitude of the maximum or minimum z of the symmeterized +table) must be greater than the full width of the bunch (excluding 1\% outlier particles). If not, +an error message is generated and tracking is stopped. + +The optional \vn{position_dependence} component determines if the wake depends upon the transverse coordinates +of the bunch particles. See \sref{s:sr.wake.syntax} for details. The default is \vn{none}. + +The wake is applied to a bunch by convoluting the bunch distribution with the wake. The optional +\vn{smoothing_sigma} component, if set non-zero, applies a Gaussian smoothing filter to the +convolution. The value of \vn{smoothing_sigma} is in meters if \vn{time_based} is False or seconds +otherwise. + +The \vn{z_scale} and \vn{amp_scale} parameters of the \vn{sr_wake} structure (\sref{s:sr.wake.syntax}) +are used with \vn{z_long}. + %----------------------------------------------------------------- \subsection{Short-Range Wakes --- Old Format} \label{s:sr.wake.old} diff --git a/bmad/doc/cover-page.tex b/bmad/doc/cover-page.tex index f275a0c02c..d56792147e 100644 --- a/bmad/doc/cover-page.tex +++ b/bmad/doc/cover-page.tex @@ -3,7 +3,7 @@ \begin{flushright} \large - Revision: November 9, 2024 \\ + Revision: November 14, 2024 \\ \end{flushright} \pdfbookmark[0]{Preamble}{Preamble} diff --git a/bmad/modules/bmad_struct.f90 b/bmad/modules/bmad_struct.f90 index c00e9b21d8..b2545bb98f 100644 --- a/bmad/modules/bmad_struct.f90 +++ b/bmad/modules/bmad_struct.f90 @@ -601,9 +601,10 @@ module bmad_struct complex(rp), allocatable :: fbunch(:), w_out(:) ! Scratch space. real(rp) :: dz = 0 ! Distance between points. real(rp) :: z0 = 0 ! Wake extent is [-z0, z0]. - integer :: plane = not_set$ ! x$, y$, xy$, z$. - integer :: position_dependence = not_set$ ! Transverse: leading$, trailing$, none$ + real(rp) :: smoothing_sigma = 0 ! 0 => No smoothing. + integer :: position_dependence = none$ ! Transverse: leading$, trailing$, none$ ! Longitudinal: x_leading$, ..., y_trailing$, none$ + logical :: time_based = .false. ! Was input time based? end type type wake_sr_mode_struct ! Psudo-mode Short-range wake struct diff --git a/bmad/modules/equality_mod.f90 b/bmad/modules/equality_mod.f90 index ea9171ac1d..b636f104ad 100644 --- a/bmad/modules/equality_mod.f90 +++ b/bmad/modules/equality_mod.f90 @@ -448,10 +448,12 @@ elemental function eq_wake_sr_z_long (f1, f2) result (is_eq) is_eq = is_eq .and. (f1%dz == f2%dz) !! f_side.equality_test[real, 0, NOT] is_eq = is_eq .and. (f1%z0 == f2%z0) -!! f_side.equality_test[integer, 0, NOT] -is_eq = is_eq .and. (f1%plane == f2%plane) +!! f_side.equality_test[real, 0, NOT] +is_eq = is_eq .and. (f1%smoothing_sigma == f2%smoothing_sigma) !! f_side.equality_test[integer, 0, NOT] is_eq = is_eq .and. (f1%position_dependence == f2%position_dependence) +!! f_side.equality_test[logical, 0, NOT] +is_eq = is_eq .and. (f1%time_based .eqv. f2%time_based) end function eq_wake_sr_z_long diff --git a/bmad/multiparticle/wake_mod.f90 b/bmad/multiparticle/wake_mod.f90 index 89bfca1957..142755b43a 100644 --- a/bmad/multiparticle/wake_mod.f90 +++ b/bmad/multiparticle/wake_mod.f90 @@ -482,22 +482,21 @@ end subroutine sr_transverse_wake_particle !-------------------------------------------------------------------------- !-------------------------------------------------------------------------- !+ -! Subroutine sr_z_long_wake (ele, bunch, z_min, z_max) +! Subroutine sr_z_long_wake (ele, bunch, z_ave) ! ! Subroutine to apply the short-range z-wake kick to a particle. ! ! Input: ! ele -- ele_struct: Element with wake. ! bunch -- bunch_struct: Bunch before wake applied. -! z_min -- real(rp): Minimum z-position of all live particles. -! z_max -- real(rp): Maximum z-position of all live particles. +! z_ave -- real(rp): Average z-position of all live particles. ! ! Output: ! orbit -- coord_struct: Ending particle coords. ! bunch -- bunch_struct: Bunch before wake applied. !+ -subroutine sr_z_long_wake (ele, bunch, z_min, z_max) +subroutine sr_z_long_wake (ele, bunch, z_ave) use spline_mod @@ -507,10 +506,10 @@ subroutine sr_z_long_wake (ele, bunch, z_min, z_max) type (wake_sr_z_long_struct), pointer :: srz type (coord_struct), pointer :: p, orbit -real(rp) x, f0, ff, f_add, kick, dz -real(rp) z_min, z_max, z_ave +real(rp) x, f0, ff, f_add, kick, dz, rz_rel, r1, r2 +real(rp) z_ave -integer i, j, ix, n2 +integer i, j, ix1, ix2, n2, n_bad, nn logical ok character(*), parameter :: r_name = 'sr_z_long_wake' @@ -523,37 +522,51 @@ subroutine sr_z_long_wake (ele, bunch, z_min, z_max) srz => sr%z_long if (srz%dz == 0) return -if (sr%z_scale * (z_max - z_min) > srz%z0) then - call out_io (s_error$, r_name, & - 'The bunch is longer than the sr-z wake can handle for element: ' // ele%name) - p%state = lost$ - return -endif - f0 = sr%amp_scale * bunch%charge_live if (sr%scale_with_length) f0 = f0 * ele%value(l$) -! Compute wake +! Compute binned bunch distribution and wake -z_ave = 0.5_rp * (z_min + z_max) -n2 = size(srz%w_out) / 2 +nn = size(srz%w_out) +n2 = (nn - 1) / 2 srz%w_out = 0 +n_bad = 0 do i = 1, size(bunch%particle) p => bunch%particle(i) - ix = nint(sr%z_scale * (p%vec(5) - z_ave) / srz%dz) + n2 if (p%state /= alive$) cycle - select case (srz%plane) + rz_rel = sr%z_scale * (p%vec(5) - z_ave) / srz%dz + n2 + 1 + ix1 = floor(rz_rel) + ix2 = ceiling(rz_rel) + if (ix1 < 1 .or. ix2 > nn) then + n_bad = n_bad + 1 + cycle + endif + + r1 = ix2 - rz_rel + r2 = rz_rel - ix1 + + select case (srz%position_dependence) case (none$, x_trailing$, y_trailing$) - srz%w_out(ix) = srz%w_out(ix) + 1 + srz%w_out(ix1) = srz%w_out(ix1) + r1 + srz%w_out(ix2) = srz%w_out(ix2) + r2 case (x_leading$) - srz%w_out(ix) = srz%w_out(ix) + p%vec(1) + srz%w_out(ix1) = srz%w_out(ix1) + r1 * p%vec(1) + srz%w_out(ix2) = srz%w_out(ix2) + r2 * p%vec(1) case (y_leading$) - srz%w_out(ix) = srz%w_out(ix) + p%vec(3) + srz%w_out(ix1) = srz%w_out(ix2) + r1 * p%vec(3) + srz%w_out(ix2) = srz%w_out(ix2) + r2 * p%vec(3) end select enddo +if (n_bad > 0.01 * size(bunch%particle)) then + call out_io (s_error$, r_name, & + 'The bunch is longer than the sr-z wake can handle for element: ' // ele%name) + p%state = lost$ + return +endif + call fft_1d(srz%w_out, -1) srz%w_out = srz%w_out * srz%fw * f0 call fft_1d(srz%w_out, 1) @@ -563,36 +576,21 @@ subroutine sr_z_long_wake (ele, bunch, z_min, z_max) do i = 1, size(bunch%particle) p => bunch%particle(i) if (p%state /= alive$) cycle - ix = nint((p%vec(5) - z_ave) / srz%dz) + n2 - select case (srz%plane) - case (none$, x_leading$, y_leading$) - p%vec(6) = p%vec(6) - srz%w_out(ix) - case (x_trailing$) - p%vec(6) = p%vec(6) - srz%w_out(ix) * p%vec(1) - case (y_trailing$) - p%vec(6) = p%vec(6) - srz%w_out(ix) * p%vec(3) - end select -enddo + rz_rel = sr%z_scale * (p%vec(5) - z_ave) / srz%dz + n2 + 1 + ix1 = floor(rz_rel) + ix2 = ceiling(rz_rel) -call fft_1d(srz%w_out, -1) -srz%w_out = srz%w_out * srz%fw * f0 -call fft_1d(srz%w_out, 1) - -! Apply wake - -do i = 1, size(bunch%particle) - p => bunch%particle(i) - if (p%state /= alive$) cycle - ix = nint(sr%z_scale * (p%vec(5) - z_ave) / srz%dz) + n2 + r1 = ix2 - rz_rel + r2 = rz_rel - ix1 - select case (srz%plane) + select case (srz%position_dependence) case (none$, x_leading$, y_leading$) - p%vec(6) = p%vec(6) - srz%w_out(ix) + p%vec(6) = p%vec(6) - (r1 * srz%w_out(ix1) + r2 * srz%w_out(ix2)) case (x_trailing$) - p%vec(6) = p%vec(6) - srz%w_out(ix) * p%vec(1) + p%vec(6) = p%vec(6) - (r1 * srz%w_out(ix1) + r2 * srz%w_out(ix2)) * p%vec(1) case (y_trailing$) - p%vec(6) = p%vec(6) - srz%w_out(ix) * p%vec(3) + p%vec(6) = p%vec(6) - (r1 * srz%w_out(ix1) + r2 * srz%w_out(ix2)) * p%vec(3) end select enddo @@ -751,7 +749,7 @@ subroutine track1_sr_wake (bunch, ele) ! Z-wake -call sr_z_long_wake(ele, bunch, p(i1)%vec(5), p(i2)%vec(5)) +call sr_z_long_wake(ele, bunch, p((i1+i2)/2)%vec(5)) ! Loop over all particles in the bunch and apply the mode wakes diff --git a/bmad/output/type_ele.f90 b/bmad/output/type_ele.f90 index 5185e6eb21..77b758b998 100644 --- a/bmad/output/type_ele.f90 +++ b/bmad/output/type_ele.f90 @@ -1422,14 +1422,10 @@ subroutine type_ele (ele, type_zero_attrib, type_mat6, type_taylor, twiss_out, t call re_allocate (li, nl+size(ele%wake%sr%z_long%w)+100, .false.) nl=nl+1; li(nl) = ' Short-Range Z-dependent Longitudinal wake:' srz => ele%wake%sr%z_long - if (srz%plane == z$) then - nl=nl+1; li(nl) = ' plane = ' // trim(sr_z_plane_name(srz%plane)) - else - nl=nl+1; li(nl) = ' plane = ' // trim(sr_z_plane_name(srz%plane)) // & - ', position_dependence = ' // trim(sr_transverse_position_dep_name(srz%position_dependence)) - endif - nl=nl+1; li(nl) = ' dz = ' // to_str(srz%dz) - nl=nl+1; li(nl) = ' +/- Wake range: ' // to_str(srz%z0) + nl=nl+1; li(nl) = ' smoothing_sigma = ' // to_str(srz%smoothing_sigma) + nl=nl+1; li(nl) = ' position_dependence = ' // trim(sr_transverse_position_dep_name(srz%position_dependence)) + nl=nl+1; li(nl) = ' dz = ' // to_str(srz%dz) + nl=nl+1; li(nl) = ' +/- Wake range: ' // to_str(srz%z0) nl=nl+1; write (li(nl), '(a, i0)') ' # wake points: ', size(srz%w) else diff --git a/bmad/parsing/bmad_parser_mod.f90 b/bmad/parsing/bmad_parser_mod.f90 index 5a24b6cfc9..7c0f39cd5e 100644 --- a/bmad/parsing/bmad_parser_mod.f90 +++ b/bmad/parsing/bmad_parser_mod.f90 @@ -1371,10 +1371,11 @@ subroutine parser_read_sr_wake (ele, delim, delim_found, err_flag) type (wake_sr_z_long_struct), pointer :: srz type (wake_sr_struct), pointer :: wake_sr +real(rp) f real(rp), allocatable :: table(:,:) integer i, itrans, ilong, iz, ipt, ix_word, n0, n1, nn, nt -logical delim_found, err_flag, err, time_based +logical delim_found, err_flag, err character(40) err_str, word, attrib_name character(1) delim @@ -1390,8 +1391,8 @@ subroutine parser_read_sr_wake (ele, delim, delim_found, err_flag) trans = wake_sr_mode_struct() long = wake_sr_mode_struct() err_flag = .true. -time_based = .false. srz => wake_sr%z_long +srz%time_based = .false. ! get data @@ -1458,13 +1459,13 @@ subroutine parser_read_sr_wake (ele, delim, delim_found, err_flag) if (.not. expect_one_of (',}', .false., ele%name, delim, delim_found)) return - case ('PLANE') - call get_switch ('SR_WAKE Z PLANE', sr_z_plane_name, srz%plane, err, ele, delim, delim_found); if (err) return + case ('SMOOTHING_SIGMA') + call parse_evaluate_value (err_str, srz%smoothing_sigma, lat, delim, delim_found, err_flag, ',', ele); if (err_flag) return case ('POSITION_DEPENDENCE') call get_switch ('SR_WAKE Z POSITION_DEPENDENCE', sr_longitudinal_position_dep_name, srz%position_dependence, err_flag, ele, delim, delim_found) if (err_flag) return - case ('T_BASED') - call parser_get_logical (attrib_name, time_based, ele%name, delim, delim_found, err_flag); if (err_flag) return + case ('TIME_BASED') + call parser_get_logical (attrib_name, srz%time_based, ele%name, delim, delim_found, err_flag); if (err_flag) return case default call parser_error ('UNKNOWN SR_WAKE Z COMPONENT: ' // attrib_name, 'FOR ELEMENT: ' // ele%name) return @@ -1524,12 +1525,24 @@ subroutine parser_read_sr_wake (ele, delim, delim_found, err_flag) allocate (wake_sr%trans(itrans)) wake_sr%trans = trans(1:itrans) -if (.not. allocated(srz%w)) then - allocate (srz%w(0), srz%fw(0), srz%w_out(0), srz%fbunch(0)) -else - if (time_based) srz%dz = srz%dz * c_light +if (allocated(srz%w)) then + if (srz%time_based) then + srz%dz = c_light * srz%dz + srz%z0 = c_light * srz%z0 + srz%w = c_light * srz%w(nt:1:-1) + srz%smoothing_sigma = c_light * srz%smoothing_sigma + endif srz%fw = srz%w call fft_1d(srz%fw, -1) + if (srz%smoothing_sigma /= 0) then + do i = 1, nt + f = real(i - 1, rp) / (nt - 1) + srz%fw(i) = srz%fw(i) * exp(-2*pi*(f*srz%smoothing_sigma)**2) + enddo + endif + +else + allocate (srz%w(0), srz%fw(0), srz%w_out(0), srz%fbunch(0)) endif err_flag = .false. diff --git a/bmad/parsing/read_digested_bmad_file.f90 b/bmad/parsing/read_digested_bmad_file.f90 index fd4ce37250..c93bb3be6a 100644 --- a/bmad/parsing/read_digested_bmad_file.f90 +++ b/bmad/parsing/read_digested_bmad_file.f90 @@ -495,6 +495,7 @@ subroutine read_this_ele (ele, ix_ele_in, error) type (converter_prob_pc_r_struct), pointer :: ppcr type (converter_direction_out_struct), pointer :: c_dir type (control_ramp1_struct), pointer ::rmp +type (wake_sr_z_long_struct), pointer :: srz integer i, j, lb1, lb2, lb3, ub1, ub2, ub3, n_cyl, n_cart, n_gen, n_grid, ix_ele, ix_branch, ix_wall3d integer i_min(3), i_max(3), ix_ele_in, ix_t(6), ios, k_max, ix_e, n_angle, n_energy @@ -913,9 +914,10 @@ subroutine read_this_ele (ele, ix_ele_in, error) read (d_unit, err = 9800, end = 9800) wake%sr%trans(i) enddo - read (d_unit, err = 9800, end = 9800) wake%sr%z_long%plane, wake%sr%z_long%position_dependence, wake%sr%z_long%dz, wake%sr%z_long%z0 - do i = 1, size(wake%sr%z_long%w) - read (d_unit, err = 9800, end = 9800) wake%sr%z_long%w(i), wake%sr%z_long%fw(i) + srz => wake%sr%z_long + read (d_unit, err = 9800, end = 9800) srz%smoothing_sigma, srz%position_dependence, srz%dz, srz%z0, srz%time_based + do i = 1, size(srz%w) + read (d_unit, err = 9800, end = 9800) srz%w(i), srz%fw(i) enddo read (d_unit, err = 9800, end = 9800) wake%lr%t_ref, wake%lr%freq_spread, wake%lr%self_wake_on, wake%lr%amp_scale, wake%lr%time_scale diff --git a/bmad/parsing/write_digested_bmad_file.f90 b/bmad/parsing/write_digested_bmad_file.f90 index 84eef34ee5..6c7a31d7af 100644 --- a/bmad/parsing/write_digested_bmad_file.f90 +++ b/bmad/parsing/write_digested_bmad_file.f90 @@ -228,6 +228,7 @@ subroutine write_this_ele (ele) type (converter_prob_pc_r_struct), pointer :: ppcr type (converter_direction_out_struct), pointer :: c_dir type (control_ramp1_struct), pointer ::rmp +type (wake_sr_z_long_struct), pointer :: srz integer ix_wall3d, ix_r, ix_d, ix_m, ix_e, ix_t(6), ix_st(0:3), ie, ib, ix_wall3d_branch integer ix_sr_long, ix_sr_trans, ix_sr_z, ix_lr_mode, ie_max, ix_s, n_var, ix_ptr, im, n1, n2 @@ -660,9 +661,10 @@ subroutine write_this_ele (ele) write (d_unit) wake%sr%trans(i) enddo - write (d_unit) wake%sr%z_long%plane, wake%sr%z_long%position_dependence, wake%sr%z_long%dz, wake%sr%z_long%z0 - do i = 1, size(wake%sr%z_long%w) - write (d_unit) wake%sr%z_long%w(i), wake%sr%z_long%fw(i) + srz => wake%sr%z_long + write (d_unit) srz%smoothing_sigma, srz%position_dependence, srz%dz, srz%z0, srz%time_based + do i = 1, size(srz%w) + write (d_unit) srz%w(i), srz%fw(i) enddo write (d_unit) wake%lr%t_ref, wake%lr%freq_spread, wake%lr%self_wake_on, wake%lr%amp_scale, wake%lr%time_scale diff --git a/cpp_bmad_interface/code/bmad_cpp_convert_mod.f90 b/cpp_bmad_interface/code/bmad_cpp_convert_mod.f90 index 995a366898..c9d0515dae 100644 --- a/cpp_bmad_interface/code/bmad_cpp_convert_mod.f90 +++ b/cpp_bmad_interface/code/bmad_cpp_convert_mod.f90 @@ -2199,14 +2199,15 @@ subroutine wake_sr_z_long_to_c (Fp, C) bind(c) interface !! f_side.to_c2_f2_sub_arg subroutine wake_sr_z_long_to_c2 (C, z_w, n1_w, z_fw, n1_fw, z_fbunch, n1_fbunch, z_w_out, & - n1_w_out, z_dz, z_z0, z_plane, z_position_dependence) bind(c) + n1_w_out, z_dz, z_z0, z_smoothing_sigma, z_position_dependence, z_time_based) bind(c) import c_bool, c_double, c_ptr, c_char, c_int, c_long, c_double_complex !! f_side.to_c2_type :: f_side.to_c2_name type(c_ptr), value :: C - real(c_double) :: z_w(*), z_dz, z_z0 + real(c_double) :: z_w(*), z_dz, z_z0, z_smoothing_sigma integer(c_int), value :: n1_w, n1_fw, n1_fbunch, n1_w_out complex(c_double_complex) :: z_fw(*), z_fbunch(*), z_w_out(*) - integer(c_int) :: z_plane, z_position_dependence + integer(c_int) :: z_position_dependence + logical(c_bool) :: z_time_based end subroutine end interface @@ -2248,7 +2249,7 @@ subroutine wake_sr_z_long_to_c2 (C, z_w, n1_w, z_fw, n1_fw, z_fbunch, n1_fbunch, !! f_side.to_c2_call call wake_sr_z_long_to_c2 (C, fvec2vec(F%w, n1_w), n1_w, fvec2vec(F%fw, n1_fw), n1_fw, & fvec2vec(F%fbunch, n1_fbunch), n1_fbunch, fvec2vec(F%w_out, n1_w_out), n1_w_out, F%dz, & - F%z0, F%plane, F%position_dependence) + F%z0, F%smoothing_sigma, F%position_dependence, c_logic(F%time_based)) end subroutine wake_sr_z_long_to_c @@ -2269,7 +2270,7 @@ end subroutine wake_sr_z_long_to_c !! f_side.to_c2_f2_sub_arg subroutine wake_sr_z_long_to_f2 (Fp, z_w, n1_w, z_fw, n1_fw, z_fbunch, n1_fbunch, z_w_out, & - n1_w_out, z_dz, z_z0, z_plane, z_position_dependence) bind(c) + n1_w_out, z_dz, z_z0, z_smoothing_sigma, z_position_dependence, z_time_based) bind(c) implicit none @@ -2282,8 +2283,9 @@ subroutine wake_sr_z_long_to_f2 (Fp, z_w, n1_w, z_fw, n1_fw, z_fbunch, n1_fbunch real(c_double), pointer :: f_w(:) integer(c_int), value :: n1_w, n1_fw, n1_fbunch, n1_w_out complex(c_double_complex), pointer :: f_fw(:), f_fbunch(:), f_w_out(:) -real(c_double) :: z_dz, z_z0 -integer(c_int) :: z_plane, z_position_dependence +real(c_double) :: z_dz, z_z0, z_smoothing_sigma +integer(c_int) :: z_position_dependence +logical(c_bool) :: z_time_based call c_f_pointer (Fp, F) @@ -2343,10 +2345,12 @@ subroutine wake_sr_z_long_to_f2 (Fp, z_w, n1_w, z_fw, n1_fw, z_fbunch, n1_fbunch F%dz = z_dz !! f_side.to_f2_trans[real, 0, NOT] F%z0 = z_z0 -!! f_side.to_f2_trans[integer, 0, NOT] -F%plane = z_plane +!! f_side.to_f2_trans[real, 0, NOT] +F%smoothing_sigma = z_smoothing_sigma !! f_side.to_f2_trans[integer, 0, NOT] F%position_dependence = z_position_dependence +!! f_side.to_f2_trans[logical, 0, NOT] +F%time_based = f_logic(z_time_based) end subroutine wake_sr_z_long_to_f2 diff --git a/cpp_bmad_interface/code/cpp_bmad_convert.cpp b/cpp_bmad_interface/code/cpp_bmad_convert.cpp index 085fcff973..2061e30c36 100644 --- a/cpp_bmad_interface/code/cpp_bmad_convert.cpp +++ b/cpp_bmad_interface/code/cpp_bmad_convert.cpp @@ -553,7 +553,8 @@ extern "C" void wake_sr_z_long_to_c (const Opaque_wake_sr_z_long_class*, CPP_wak // c_side.to_f2_arg extern "C" void wake_sr_z_long_to_f2 (Opaque_wake_sr_z_long_class*, c_RealArr, Int, - c_ComplexArr, Int, c_ComplexArr, Int, c_ComplexArr, Int, c_Real&, c_Real&, c_Int&, c_Int&); + c_ComplexArr, Int, c_ComplexArr, Int, c_ComplexArr, Int, c_Real&, c_Real&, c_Real&, c_Int&, + c_Bool&); extern "C" void wake_sr_z_long_to_f (const CPP_wake_sr_z_long& C, Opaque_wake_sr_z_long_class* F) { // c_side.to_f_setup[real, 1, ALLOC] @@ -583,14 +584,15 @@ extern "C" void wake_sr_z_long_to_f (const CPP_wake_sr_z_long& C, Opaque_wake_sr // c_side.to_f2_call wake_sr_z_long_to_f2 (F, z_w, n1_w, z_fw, n1_fw, z_fbunch, n1_fbunch, z_w_out, n1_w_out, - C.dz, C.z0, C.plane, C.position_dependence); + C.dz, C.z0, C.smoothing_sigma, C.position_dependence, C.time_based); } // c_side.to_c2_arg extern "C" void wake_sr_z_long_to_c2 (CPP_wake_sr_z_long& C, c_RealArr z_w, Int n1_w, c_ComplexArr z_fw, Int n1_fw, c_ComplexArr z_fbunch, Int n1_fbunch, c_ComplexArr z_w_out, - Int n1_w_out, c_Real& z_dz, c_Real& z_z0, c_Int& z_plane, c_Int& z_position_dependence) { + Int n1_w_out, c_Real& z_dz, c_Real& z_z0, c_Real& z_smoothing_sigma, c_Int& + z_position_dependence, c_Bool& z_time_based) { // c_side.to_c2_set[real, 1, ALLOC] @@ -616,10 +618,12 @@ extern "C" void wake_sr_z_long_to_c2 (CPP_wake_sr_z_long& C, c_RealArr z_w, Int C.dz = z_dz; // c_side.to_c2_set[real, 0, NOT] C.z0 = z_z0; - // c_side.to_c2_set[integer, 0, NOT] - C.plane = z_plane; + // c_side.to_c2_set[real, 0, NOT] + C.smoothing_sigma = z_smoothing_sigma; // c_side.to_c2_set[integer, 0, NOT] C.position_dependence = z_position_dependence; + // c_side.to_c2_set[logical, 0, NOT] + C.time_based = z_time_based; } //-------------------------------------------------------------------- diff --git a/cpp_bmad_interface/code/cpp_equality.cpp b/cpp_bmad_interface/code/cpp_equality.cpp index a17446d074..8ecb35e4e4 100644 --- a/cpp_bmad_interface/code/cpp_equality.cpp +++ b/cpp_bmad_interface/code/cpp_equality.cpp @@ -269,8 +269,9 @@ bool operator== (const CPP_wake_sr_z_long& x, const CPP_wake_sr_z_long& y) { is_eq = is_eq && is_all_equal(x.w_out, y.w_out); is_eq = is_eq && (x.dz == y.dz); is_eq = is_eq && (x.z0 == y.z0); - is_eq = is_eq && (x.plane == y.plane); + is_eq = is_eq && (x.smoothing_sigma == y.smoothing_sigma); is_eq = is_eq && (x.position_dependence == y.position_dependence); + is_eq = is_eq && (x.time_based == y.time_based); return is_eq; }; diff --git a/cpp_bmad_interface/include/cpp_bmad_classes.h b/cpp_bmad_interface/include/cpp_bmad_classes.h index 344313c23b..e5a1c6e770 100644 --- a/cpp_bmad_interface/include/cpp_bmad_classes.h +++ b/cpp_bmad_interface/include/cpp_bmad_classes.h @@ -896,8 +896,9 @@ class CPP_wake_sr_z_long { Complex_ARRAY w_out; Real dz; Real z0; - Int plane; + Real smoothing_sigma; Int position_dependence; + Bool time_based; CPP_wake_sr_z_long() : w(0.0, 0), @@ -906,8 +907,9 @@ class CPP_wake_sr_z_long { w_out(0.0, 0), dz(0.0), z0(0.0), - plane(Bmad::NOT_SET), - position_dependence(Bmad::NOT_SET) + smoothing_sigma(0.0), + position_dependence(Bmad::NONE), + time_based(false) {} ~CPP_wake_sr_z_long() { diff --git a/cpp_bmad_interface/interface_test/bmad_cpp_test_mod.f90 b/cpp_bmad_interface/interface_test/bmad_cpp_test_mod.f90 index d367ae8898..f36ff38ddc 100644 --- a/cpp_bmad_interface/interface_test/bmad_cpp_test_mod.f90 +++ b/cpp_bmad_interface/interface_test/bmad_cpp_test_mod.f90 @@ -1392,10 +1392,12 @@ subroutine set_wake_sr_z_long_test_pattern (F, ix_patt) rhs = 9 + offset; F%dz = rhs !! f_side.test_pat[real, 0, NOT] rhs = 10 + offset; F%z0 = rhs -!! f_side.test_pat[integer, 0, NOT] -rhs = 11 + offset; F%plane = rhs +!! f_side.test_pat[real, 0, NOT] +rhs = 11 + offset; F%smoothing_sigma = rhs !! f_side.test_pat[integer, 0, NOT] rhs = 12 + offset; F%position_dependence = rhs +!! f_side.test_pat[logical, 0, NOT] +rhs = 13 + offset; F%time_based = (modulo(rhs, 2) == 0) end subroutine set_wake_sr_z_long_test_pattern diff --git a/cpp_bmad_interface/interface_test/cpp_bmad_test.cpp b/cpp_bmad_interface/interface_test/cpp_bmad_test.cpp index 1aba32d065..58a881c8f2 100644 --- a/cpp_bmad_interface/interface_test/cpp_bmad_test.cpp +++ b/cpp_bmad_interface/interface_test/cpp_bmad_test.cpp @@ -900,12 +900,15 @@ void set_CPP_wake_sr_z_long_test_pattern (CPP_wake_sr_z_long& C, int ix_patt) { // c_side.test_pat[real, 0, NOT] rhs = 10 + offset; C.z0 = rhs; - // c_side.test_pat[integer, 0, NOT] - rhs = 11 + offset; C.plane = rhs; + // c_side.test_pat[real, 0, NOT] + rhs = 11 + offset; C.smoothing_sigma = rhs; // c_side.test_pat[integer, 0, NOT] rhs = 12 + offset; C.position_dependence = rhs; + // c_side.test_pat[logical, 0, NOT] + rhs = 13 + offset; C.time_based = (rhs % 2 == 0); + } diff --git a/tao/version/tao_version_mod.f90 b/tao/version/tao_version_mod.f90 index b73cd827d0..f831189552 100644 --- a/tao/version/tao_version_mod.f90 +++ b/tao/version/tao_version_mod.f90 @@ -6,5 +6,5 @@ !- module tao_version_mod -character(*), parameter :: tao_version_date = "2024/11/10 17:17:09" +character(*), parameter :: tao_version_date = "2024/11/12 01:38:51" end module diff --git a/util_programs/wake_plot/README b/util_programs/wake_plot/README index d6f37e7de3..55a264d85f 100644 --- a/util_programs/wake_plot/README +++ b/util_programs/wake_plot/README @@ -51,7 +51,7 @@ An input file uses Fortan namelist syntax. Parameters are read from a namelist n which looks like: ¶ms - lat_file = "lat.bmad" + lat_file = "lat.bmad" ... etc... / @@ -67,18 +67,18 @@ lat_file: Bmad lattice file containing the element with the wakes. wake_ele_name: Wake element identifier. Can be name or element index. Must match to a unique element. who: What wake type to use. Possibilities are: - "lr" -- Long range wakes. - "sr" -- Short range wakes. + "lr" -- Long range wakes. + "sr" -- Short range wakes. "sr-mode" -- Short range pseudo-mode wakes. Z-position wakes are ignored. "sr-long" -- Short range longitudinal pseudo-mode wake. Only used with ix_wake non-zero. "sr-trans" -- Short range transverse pseudo-mode wake. Only used with ix_wake non-zero. - "sr-z" -- Short range Z-position wakes. Pseudo-mode wakes are ignored. + "sr-z-long" -- Short range Z-position wakes. Pseudo-mode wakes are ignored. When multiple modes are present, the sum of the wakes is plotted. m_order: Mode order to plot. Only used with LR wakes. m_order = 1 => dipole wakes, etc. Ignored if ix_wake is nonzero. Must be set if there is more than one wake and ix_wake = 0. -ix_wake: Wake index to plot. Only used with who set to "lr", "sr-long", "sr-trans", or "sr-z" where +ix_wake: Wake index to plot. Only used with "who" set to "lr", "sr-long", or "sr-trans" where there is a single list of wakes. plot only the wake with the index corresponding to ix_wake (ix_wake = 1 --> Only use first wake in list, etc). Default is zero which means ignore ix_order. If non-zero, overrides any setting of m_order. @@ -99,13 +99,14 @@ make_plot: If False, do not make a plot (used if only a postscript file is wante postscript_file: If non-blank, Make a postscript file. -Plot_limit: For LR plots, the wake is plotted over a time range [0, plot_limit]. For SR plots, -the wake is plotted over a z range [plot_limit, 0] where plot_limit should be negative in this case. +Plot_limit: + * For LR plots, the wake is plotted over a time range [0, plot_limit]. + * For "sr-z-long" plot the z range is [-plot_limit, +plot_limit] except if plot_limit is zero then the entire + wake table is plotted. + * For SR plots, the wake is plotted over a z range [plot_limit, 0] where plot_limit should be negative in this case. -n_points: Number of data points to use in drawing the wake. - -draw_knot_points: Only used with SR wakes. If True, draw symbols on the wake curves at the z-positions -corresponding to the spline knot points of the Z-position wakes. +n_points: Number of data points to use in drawing the wake. Not used with "sr-z-long" where the number +of points is determined by the number points in the wake table that are visible in the plot text_scale: Scale for drawing text. Default is 1.2 diff --git a/util_programs/wake_plot/example/wake_plot.init b/util_programs/wake_plot/example/wake_plot.init index 7295702dfe..861e30b835 100644 --- a/util_programs/wake_plot/example/wake_plot.init +++ b/util_programs/wake_plot/example/wake_plot.init @@ -1,6 +1,6 @@ ¶ms lat_file = "wake.bmad" ! Bmad lattice file. - wake_ele_name = "p1" ! Lattice element with wakes. + wake_ele_name = "p1" ! Name or index of lattice element with wakes. who = "sr" ! What wake type to use. m_order = 0 ! Order of LR wake to use. 0 => all orders. ix_wake = 0 ! Index of wake to use. 0 => Ignore this parameter. diff --git a/util_programs/wake_plot/wake_plot.f90 b/util_programs/wake_plot/wake_plot.f90 index 67b9ae1339..45a9a43529 100644 --- a/util_programs/wake_plot/wake_plot.f90 +++ b/util_programs/wake_plot/wake_plot.f90 @@ -11,15 +11,16 @@ program wake_plot type (ele_pointer_struct), allocatable, target :: eles(:) type (wake_sr_struct), pointer :: sr type (wake_lr_struct), pointer :: lr +type (wake_sr_z_long_struct), pointer :: srz type (wake_lr_mode_struct), allocatable :: lr_mode(:) type (bunch_struct), target :: bunch type (coord_struct), pointer :: p1, p2 real(rp) :: plot_min, plot_max, plot_limit, plot_size(2), text_scale, zero6(6) = 0, vec2(6) real(rp) xy_leading(2), xy_trailing(2), leading_charge, z -real(rp), allocatable :: time(:), z_z(:), w_z(:), wx(:), wy(:), wz(:) +real(rp), allocatable :: x_axis(:), wx(:), wy(:), wz(:) -integer n, n_points, ix_wake, m_order, i, im, iz, ik, n_loc +integer n, n_points, ix_wake, m_order, i, im, iz, ik, n_loc, ix, nw, n0 logical make_plot, err, err_flag, ok @@ -49,8 +50,6 @@ program wake_plot read (1, nml = params) close (1) -allocate (time(n_points), wx(n_points), wy(n_points), wz(n_points)) - ! Init call bmad_parser(lat_file, lat, err_flag = err) @@ -94,10 +93,41 @@ program wake_plot p1%charge = leading_charge endif +! + +if (who == 'sr-z-long') then + srz => wake_ele%wake%sr%z_long + if (size(srz%w) == 0) then + print *, 'No SR Z_long wake.' + stop + endif + + if (plot_limit == 0 .or. plot_limit > srz%z0) then + plot_limit= srz%z0 + endif + + nw = (size(srz%w) - 1) / 2 + n0 = nint(plot_limit / srz%dz) + allocate(wz(2*n0+1), x_axis(2*n0+1)) + plot_min = -plot_limit + plot_max = plot_limit + + do im = -n0, n0 + ix = im + n0 + 1 + x_axis(ix) = im * srz%dz + wz(ix) = srz%w(im + nw + 1) + enddo + + if (make_plot) call make_this_plot ('X') + if (postscript_file /= '') call make_this_plot ('PS') + stop +endif + ! -select case (who(1:2)) +allocate (x_axis(n_points), wx(n_points), wy(n_points), wz(n_points)) +select case (who(1:2)) case ('lr') lr => wake_ele%wake%lr bmad_com%lr_wakes_on = .true. @@ -131,7 +161,7 @@ program wake_plot z = plot_min + (plot_max - plot_min) * (i - 1.0_rp) / (n_points - 1.0_rp) p2%vec = vec2 p2%vec(5) = z - time(i) = z + x_axis(i) = z call zero_lr_wakes_in_lat(lat) call track1_lr_wake(bunch, wake_ele) wx(i) = p2%vec(2) @@ -158,10 +188,9 @@ program wake_plot allocate(sr%z_long%w(0)) allocate(sr%long(0)) if (ix_wake /= 0) sr%trans = [sr%trans(ix_wake)] - case ('sr-z') + case ('sr-z-long') allocate(sr%trans(0)) allocate(sr%long(0)) - if (ix_wake /= 0) sr%z_long%w = [sr%z_long%w(ix_wake)] end select p2%vec(5) = -1 @@ -171,15 +200,13 @@ program wake_plot z = plot_min + (plot_max - plot_min) * (i - 1.0_rp) / (n_points - 1.0_rp) p2%vec = vec2 p2%vec(5) = z - time(i) = z + x_axis(i) = z call track1_sr_wake(bunch, wake_ele) wx(i) = p2%vec(2) wy(i) = p2%vec(4) wz(i) = p2%vec(6) enddo - allocate(z_z(0), w_z(0)) - end select ! Plotting @@ -190,19 +217,19 @@ program wake_plot !------------------------------------------------------ contains -subroutine make_this_plot(who) +subroutine make_this_plot(plot_window) type (qp_line_struct), allocatable :: line(:) type (qp_symbol_struct), allocatable :: symbol(:) real(rp) height -character(*) who +character(*) plot_window character(16) :: what_str, x_axis_str integer j, k, n, id ! -if (who == 'X') then +if (plot_window == 'X') then call qp_open_page ('X', id, plot_size(1), plot_size(2), 'POINTS') elseif (plot_size(1) > plot_size(2)) then call qp_open_page ('PS-L', id, plot_size(1), plot_size(2), 'POINTS', plot_file = postscript_file) @@ -213,9 +240,10 @@ subroutine make_this_plot(who) select case (plot_type) case ('wake'); what_str = 'Wake' case ('kick'); what_str = 'Kick' +case default; what_str = 'Wake' end select -select case (who) +select case (plot_type) case ('lr'); x_axis_str = 'Time' case default; x_axis_str = 'Z' end select @@ -227,35 +255,48 @@ subroutine make_this_plot(who) call qp_set_text_attrib ('AXIS_LABEL', text_scale*height) call qp_set_page_border (0.01_rp, 0.01_rp, 0.01_rp, 0.01_rp, '%PAGE') -call qp_set_margin (0.10_rp, 0.02_rp, 0.07_rp, 0.01_rp, '%PAGE') +call qp_set_margin (0.15_rp, 0.02_rp, 0.07_rp, 0.01_rp, '%PAGE') + +! -call qp_calc_and_set_axis ('X', plot_min, plot_max, 4, 8, 'ZERO_AT_END') +if (who == 'sr-z-long') then + call qp_calc_and_set_axis ('X', plot_min, plot_max, 4, 8, 'GENERAL') + call qp_set_box (1, 1, 1, 1) + call qp_calc_and_set_axis ('Y', minval(wz), maxval(wz), 4, 8, 'GENERAL') + call qp_draw_axes (x_axis_str, what_str) + + call qp_draw_data (x_axis, wz, .true., 0) ! -call qp_set_box (1, 3, 1, 3) -call qp_calc_and_set_axis ('Y', minval(wx), maxval(wx), 4, 8, 'GENERAL') -call qp_draw_axes (x_axis_str, what_str) +else + call qp_calc_and_set_axis ('X', plot_min, plot_max, 4, 8, 'ZERO_AT_END') -call qp_draw_data (time, wx, .true., 0) + call qp_set_box (1, 3, 1, 3) + call qp_calc_and_set_axis ('Y', minval(wx), maxval(wx), 4, 8, 'GENERAL') + call qp_draw_axes (x_axis_str, what_str) -! + call qp_draw_data (x_axis, wx, .true., 0) -call qp_set_box (1, 2, 1, 3) -call qp_calc_and_set_axis ('Y', minval(wy), maxval(wy), 4, 8, 'GENERAL') -call qp_draw_axes (x_axis_str, what_str) + ! -! + call qp_set_box (1, 2, 1, 3) + call qp_calc_and_set_axis ('Y', minval(wy), maxval(wy), 4, 8, 'GENERAL') + call qp_draw_axes (x_axis_str, what_str) + + call qp_draw_data (x_axis, wy, .true., 0) -call qp_set_box (1, 1, 1, 3) -call qp_calc_and_set_axis ('Y', minval(wz), maxval(wz), 4, 8, 'GENERAL') -call qp_draw_axes (x_axis_str, what_str) + ! -call qp_draw_data (z_z, w_z, .true., 0) + call qp_set_box (1, 1, 1, 3) + call qp_calc_and_set_axis ('Y', minval(wz), maxval(wz), 4, 8, 'GENERAL') + call qp_draw_axes (x_axis_str, what_str) + call qp_draw_data (x_axis, wz, .true., 0) +endif ! -if (who == 'X') then +if (plot_window == 'X') then write (*, '(a)', advance = 'NO') 'Hit any key to continue.' read (*, '(a)') ans endif