diff --git a/Makefile_test b/Makefile_test index f53c0d775..857877b70 100644 --- a/Makefile_test +++ b/Makefile_test @@ -60,6 +60,7 @@ test-dynap \ test-c6t test-c6t-2 test-c6t-3 test-c6t-4 \ test-thick-quad test-thick-quad-2 test-thick-quad-3 \ test-thick-dipole test-thick-dipole-2 test-thick-dipole-3 \ +test-thick-solenoid \ test-jacobian test-jacobian-2 test-jacobian-knobs \ test-match test-match-2 test-match-3 test-match-4 \ test-match-5 test-match-6 test-match-7 test-match-8 test-match-9\ diff --git a/src/trrun.f90 b/src/trrun.f90 index 59250d6e0..7aa153ddc 100644 --- a/src/trrun.f90 +++ b/src/trrun.f90 @@ -469,8 +469,12 @@ subroutine trrun(switch, turns, orbit0, rt, part_id, last_turn, last_pos, & l_buf(nlm+1) = el call element_name(el_name,len(el_name)) - if ( code.ne.code_drift .and. code.ne.code_quadrupole .and. code.ne.code_sbend & - .and. code.ne.code_matrix .and. el.ne.zero ) then ! rbend missing ? + if ( code.ne.code_drift .and. & + code.ne.code_quadrupole .and. & + code.ne.code_sbend .and. & + code.ne.code_matrix .and. & + code.ne.code_solenoid .and. & + el.ne.zero ) then ! rbend missing ? !if (.not. (is_drift() .or. is_thin() .or. is_quad() .or. is_dipole() .or. is_matrix()) ) then print *," " print *,"code: ",code," el: ",el," THICK ELEMENT FOUND" @@ -3145,7 +3149,8 @@ subroutine ttdpdg(track, ktrack) end subroutine ttdpdg subroutine trsol(track,ktrack) - use math_constfi, only : one, two, half + + use math_constfi, only : zero, one, two, half, four implicit none !----------------------------------------------------------------------* ! Purpose: * @@ -3160,66 +3165,97 @@ subroutine trsol(track,ktrack) integer :: ktrack integer :: i - double precision :: beta - double precision :: sk, skl, sks, sksl, cosTh, sinTh, Q, R, Z + double precision :: bet0, bet + double precision :: sk, skl, cosTh, sinTh, Q, R, Z double precision :: xf, yf, pxf, pyf, sigf, psigf, bvk double precision :: onedp, fpsig, fppsig double precision :: get_value, node_value + double precision :: omega, length, length_ + double precision :: x_, y_, z_, px_, py_, pt_ + !---- Initialize. - ! dtbyds = get_value('probe ','dtbyds ') - ! gamma = get_value('probe ','gamma ') - beta = get_value('probe ','beta ') - ! deltap = get_value('probe ','deltap ') - ! + bet0 = get_value('probe ','beta ') + !---- Get solenoid parameters ! elrad = node_value('lrad ') - sksl = node_value('ksi ') - sks = node_value('ks ') - - !---- BV flag bvk = node_value('other_bv ') - sks = sks * bvk - sksl = sksl * bvk - - !---- Set up strengths - ! sk = sks / two / (one + deltap) - sk = sks / two - skl = sksl / two - - !---- Loop over particles - do i = 1, ktrack - ! Ripken formulae p.28 (3.35 and 3.36) - xf = track(1,i) - yf = track(3,i) - psigf = track(6,i) / beta - - ! We do not use a constant deltap!!!!! WE use full 6D formulae! - onedp = sqrt( one + 2*psigf + (beta**2)*(psigf**2) ) - fpsig = onedp - one - fppsig = ( one + (beta**2)*psigf ) / onedp - - ! Set up C,S, Q,R,Z - cosTh = cos(skl/onedp) - sinTh = sin(skl/onedp) - Q = -skl * sk / onedp - R = fppsig / (onedp**2) * skl * sk - Z = fppsig / (onedp**2) * skl - - pxf = track(2,i) + xf*Q - pyf = track(4,i) + yf*Q - sigf = track(5,i)*beta - half*(xf**2 + yf**2)*R - - ! Ripken formulae p.29 (3.37) - track(1,i) = xf * cosTh + yf * sinTh - track(2,i) = pxf * cosTh + pyf * sinTh - track(3,i) = -xf * sinTh + yf * cosTh - track(4,i) = -pxf * sinTh + pyf * cosTh - track(5,i) = (sigf + (xf*pyf - yf*pxf)*Z) / beta - ! track(6,i) = psigf*beta - - enddo + sk = bvk * node_value('ks ') / two + length = node_value('l ') + + if (length.eq.zero) then + + skl = bvk * node_value('ksi ') / two + + !---- Loop over particles + do i = 1, ktrack + ! Ripken formulae p.28 (3.35 and 3.36) + xf = track(1,i) + yf = track(3,i) + psigf = track(6,i) / bet0 + + ! We do not use a constant deltap!!!!! WE use full 6D formulae! + onedp = sqrt( one + two*psigf + (bet0**2)*(psigf**2) ) + fpsig = onedp - one + fppsig = ( one + (bet0**2)*psigf ) / onedp + + ! Set up C,S, Q,R,Z + cosTh = cos(skl/onedp) + sinTh = sin(skl/onedp) + Q = -skl * sk / onedp + R = fppsig / (onedp**2) * skl * sk + Z = fppsig / (onedp**2) * skl + + pxf = track(2,i) + xf*Q + pyf = track(4,i) + yf*Q + sigf = track(5,i)*bet0 - half*(xf**2 + yf**2)*R + + ! Ripken formulae p.29 (3.37) + track(1,i) = xf * cosTh + yf * sinTh + track(2,i) = pxf * cosTh + pyf * sinTh + track(3,i) = -xf * sinTh + yf * cosTh + track(4,i) = -pxf * sinTh + pyf * cosTh + track(5,i) = (sigf + (xf*pyf - yf*pxf)*Z) / bet0 + ! track(6,i) = psigf*bet0 + enddo + else + if (sk.ne.zero) then + skl = sk*length + + !---- Loop over particles + do i = 1, ktrack + ! initial phase space coordinates + x_ = track(1,i) + y_ = track(3,i) + px_ = track(2,i) + py_ = track(4,i) + z_ = track(5,i) + pt_ = track(6,i) + + ! set up constants + onedp = sqrt(one + two*pt_/bet0 + pt_**2); + bet = onedp / (one/bet0 + pt_); + + ! set up constants + cosTh = cos(two*skl/onedp) + sinTh = sin(two*skl/onedp) + omega = sk/onedp; + + ! total path length traveled by the particle + length_ = length - half/(onedp**2)*(omega*(sinTh-two*length*omega)*(x_**2+y_**2)+& + two*(one-cosTh)*(px_*x_+py_*y_)-(sinTh/omega+two*length)*(px_**2+py_**2))/four; + + track(1,i) = ((one+cosTh)*x_+sinTh*y_+(px_*sinTh-py_*(cosTh-one))/omega)/two; + track(3,i) = ((one+cosTh)*y_-sinTh*x_+(py_*sinTh+px_*(cosTh-one))/omega)/two; + track(2,i) = (omega*((cosTh-one)*y_-sinTh*x_)+py_*sinTh+px_*(one+cosTh))/two; + track(4,i) = (omega*((one-cosTh)*x_-sinTh*y_)-px_*sinTh+py_*(one+cosTh))/two; + track(5,i) = z_ + length/bet0 - length_/bet; + enddo + else + call ttdrf(length,track,ktrack); + endif + endif end subroutine trsol subroutine tttrans(track,ktrack) @@ -4689,6 +4725,11 @@ subroutine trphot(el,curv,rfac,deltap) amean = five * sqrt(three) / (twelve * hbar * clight) * abs(arad * pc * (one+deltap) * el * curv) ucrit = three/two * hbar * clight * gamma**3 * abs(curv) sumxi = zero + if (amean.gt.3d-1) then + print *,"More than 0.3 photons emitted in element." + print *,"You might want to consider increasing the number of slices to reduce this number." + endif + if (amean .gt. zero) then call dpoissn(amean, nphot, ierror) diff --git a/tests/test-thick-solenoid/madx-thick.obs0001.p0001.cfg b/tests/test-thick-solenoid/madx-thick.obs0001.p0001.cfg new file mode 100644 index 000000000..a34c9f0de --- /dev/null +++ b/tests/test-thick-solenoid/madx-thick.obs0001.p0001.cfg @@ -0,0 +1,3 @@ +4-6 * skip # head + +* * any abs=1e-10 rel=5e-7 diff --git a/tests/test-thick-solenoid/madx-thick.obs0001.p0001.ref b/tests/test-thick-solenoid/madx-thick.obs0001.p0001.ref new file mode 100644 index 000000000..3ef71d4f3 --- /dev/null +++ b/tests/test-thick-solenoid/madx-thick.obs0001.p0001.ref @@ -0,0 +1,10 @@ +@ NAME %19s "TRACK.OBS0001.P0001" +@ TYPE %08s "TRACKOBS" +@ TITLE %28s "thick-solenoid tracking test" +@ ORIGIN %17s "5.03.07 Darwin 64" +@ DATE %08s "06/11/17" +@ TIME %08s "18.38.10" +* NUMBER TURN X PX Y PY T PT S E +$ %d %d %le %le %le %le %le %le %le %le + 1 0 0.001 0 0 0 0 -0.001 0 1 + 1 1 -0.0004179664169 -9.093718502e-05 0 0 -6.495195848e-08 -0.001 0 1 diff --git a/tests/test-thick-solenoid/madx-thick.obs0001.p0002.cfg b/tests/test-thick-solenoid/madx-thick.obs0001.p0002.cfg new file mode 100644 index 000000000..73364783d --- /dev/null +++ b/tests/test-thick-solenoid/madx-thick.obs0001.p0002.cfg @@ -0,0 +1,3 @@ +4-6 * skip # head + +* * any abs=1e-10 rel=5e-7 # SLC6, GCC32 diff --git a/tests/test-thick-solenoid/madx-thick.obs0001.p0002.ref b/tests/test-thick-solenoid/madx-thick.obs0001.p0002.ref new file mode 100644 index 000000000..bb9f235aa --- /dev/null +++ b/tests/test-thick-solenoid/madx-thick.obs0001.p0002.ref @@ -0,0 +1,10 @@ +@ NAME %19s "TRACK.OBS0001.P0002" +@ TYPE %08s "TRACKOBS" +@ TITLE %28s "thick-solenoid tracking test" +@ ORIGIN %17s "5.03.07 Darwin 64" +@ DATE %08s "06/11/17" +@ TIME %08s "18.38.10" +* NUMBER TURN X PX Y PY T PT S E +$ %d %d %le %le %le %le %le %le %le %le + 2 0 0.001 0 0 0 0 2.220446049e-16 0 1 + 2 1 -0.0004161468365 -9.092974268e-05 0 3.388131789e-21 -5.946003867e-08 2.220446049e-16 0 1 diff --git a/tests/test-thick-solenoid/madx-thick.obs0001.p0003.cfg b/tests/test-thick-solenoid/madx-thick.obs0001.p0003.cfg new file mode 100644 index 000000000..73364783d --- /dev/null +++ b/tests/test-thick-solenoid/madx-thick.obs0001.p0003.cfg @@ -0,0 +1,3 @@ +4-6 * skip # head + +* * any abs=1e-10 rel=5e-7 # SLC6, GCC32 diff --git a/tests/test-thick-solenoid/madx-thick.obs0001.p0003.ref b/tests/test-thick-solenoid/madx-thick.obs0001.p0003.ref new file mode 100644 index 000000000..02219aa83 --- /dev/null +++ b/tests/test-thick-solenoid/madx-thick.obs0001.p0003.ref @@ -0,0 +1,10 @@ +@ NAME %19s "TRACK.OBS0001.P0003" +@ TYPE %08s "TRACKOBS" +@ TITLE %28s "thick-solenoid tracking test" +@ ORIGIN %17s "5.03.07 Darwin 64" +@ DATE %08s "06/11/17" +@ TIME %08s "18.38.10" +* NUMBER TURN X PX Y PY T PT S E +$ %d %d %le %le %le %le %le %le %le %le + 3 0 0.001 0 0 0 0 0.001 0 1 + 3 1 -0.0004143292288 -9.092178557e-05 6.783039842e-20 0 -5.398501379e-08 0.001 0 1 diff --git a/tests/test-thick-solenoid/madx-thin.obs0001.p0001.cfg b/tests/test-thick-solenoid/madx-thin.obs0001.p0001.cfg new file mode 100644 index 000000000..a34c9f0de --- /dev/null +++ b/tests/test-thick-solenoid/madx-thin.obs0001.p0001.cfg @@ -0,0 +1,3 @@ +4-6 * skip # head + +* * any abs=1e-10 rel=5e-7 diff --git a/tests/test-thick-solenoid/madx-thin.obs0001.p0001.ref b/tests/test-thick-solenoid/madx-thin.obs0001.p0001.ref new file mode 100644 index 000000000..db9bac581 --- /dev/null +++ b/tests/test-thick-solenoid/madx-thin.obs0001.p0001.ref @@ -0,0 +1,10 @@ +@ NAME %19s "TRACK.OBS0001.P0001" +@ TYPE %08s "TRACKOBS" +@ TITLE %28s "thick-solenoid tracking test" +@ ORIGIN %17s "5.03.07 Darwin 64" +@ DATE %08s "06/11/17" +@ TIME %08s "18.38.10" +* NUMBER TURN X PX Y PY T PT S E +$ %d %d %le %le %le %le %le %le %le %le + 1 0 0.001 0 0 0 0 -0.001 0 1 + 1 1 -0.0004179968082 -9.084288303e-05 -2.711952481e-19 -3.599890026e-20 -1.054344426e-07 -0.001 0 1 diff --git a/tests/test-thick-solenoid/madx-thin.obs0001.p0002.cfg b/tests/test-thick-solenoid/madx-thin.obs0001.p0002.cfg new file mode 100644 index 000000000..73364783d --- /dev/null +++ b/tests/test-thick-solenoid/madx-thin.obs0001.p0002.cfg @@ -0,0 +1,3 @@ +4-6 * skip # head + +* * any abs=1e-10 rel=5e-7 # SLC6, GCC32 diff --git a/tests/test-thick-solenoid/madx-thin.obs0001.p0002.ref b/tests/test-thick-solenoid/madx-thin.obs0001.p0002.ref new file mode 100644 index 000000000..005978748 --- /dev/null +++ b/tests/test-thick-solenoid/madx-thin.obs0001.p0002.ref @@ -0,0 +1,10 @@ +@ NAME %19s "TRACK.OBS0001.P0002" +@ TYPE %08s "TRACKOBS" +@ TITLE %28s "thick-solenoid tracking test" +@ ORIGIN %17s "5.03.07 Darwin 64" +@ DATE %08s "06/11/17" +@ TIME %08s "18.38.10" +* NUMBER TURN X PX Y PY T PT S E +$ %d %d %le %le %le %le %le %le %le %le + 2 0 0.001 0 0 0 0 2.220446049e-16 0 1 + 2 1 -0.0004161771647 -9.092639127e-05 -2.050193423e-19 5.145725155e-20 -1.000038811e-07 2.220446049e-16 0 1 diff --git a/tests/test-thick-solenoid/madx-thin.obs0001.p0003.cfg b/tests/test-thick-solenoid/madx-thin.obs0001.p0003.cfg new file mode 100644 index 000000000..73364783d --- /dev/null +++ b/tests/test-thick-solenoid/madx-thin.obs0001.p0003.cfg @@ -0,0 +1,3 @@ +4-6 * skip # head + +* * any abs=1e-10 rel=5e-7 # SLC6, GCC32 diff --git a/tests/test-thick-solenoid/madx-thin.obs0001.p0003.ref b/tests/test-thick-solenoid/madx-thin.obs0001.p0003.ref new file mode 100644 index 000000000..bc4f106cc --- /dev/null +++ b/tests/test-thick-solenoid/madx-thin.obs0001.p0003.ref @@ -0,0 +1,10 @@ +@ NAME %19s "TRACK.OBS0001.P0003" +@ TYPE %08s "TRACKOBS" +@ TITLE %28s "thick-solenoid tracking test" +@ ORIGIN %17s "5.03.07 Darwin 64" +@ DATE %08s "06/11/17" +@ TIME %08s "18.38.10" +* NUMBER TURN X PX Y PY T PT S E +$ %d %d %le %le %le %le %le %le %le %le + 3 0 0.001 0 0 0 0 0.001 0 1 + 3 1 -0.0004143594937 -9.100936927e-05 -4.805886845e-19 -2.922263668e-20 -9.458958447e-08 0.001 0 1 diff --git a/tests/test-thick-solenoid/test-thick-solenoid.cfg b/tests/test-thick-solenoid/test-thick-solenoid.cfg new file mode 100644 index 000000000..2ae61b9be --- /dev/null +++ b/tests/test-thick-solenoid/test-thick-solenoid.cfg @@ -0,0 +1,2 @@ +1-7 * skip # head +* * any abs=1e-10 rel=5e-7 # global diff --git a/tests/test-thick-solenoid/test-thick-solenoid.madx b/tests/test-thick-solenoid/test-thick-solenoid.madx new file mode 100644 index 000000000..b77925c8b --- /dev/null +++ b/tests/test-thick-solenoid/test-thick-solenoid.madx @@ -0,0 +1,41 @@ +option, -debug, -echo; + +title, "thick-solenoid tracking test"; + +L_sol := 10.; +Ks_sol := 0.2; + +sol1: solenoid, l=L_sol, KS=Ks_sol, KSI=Ks_sol*L_sol; +sol2: solenoid, l=L_sol, KS=-Ks_sol, KSI=Ks_sol*L_sol; + +myseq: sequence, l=2*L_sol; +sol1, at:= L_sol/2; +sol2, at:= L_sol+L_sol/2; +endsequence; + +beam,energy=1; ! +use, sequence=myseq; + +track, onepass, dump, file= "madx-thick"; +start, x=1e-3, y=0, px=0, py=0, t=0, pt=-1e-3; +start, x=1e-3, y=0, px=0, py=0, t=0, pt=0.0; +start, x=1e-3, y=0, px=0, py=0, t=0, pt=1e-3; +run, turns=1; +endtrack; + +use, sequence=myseq; +select, flag=makethin, class=solenoid, thick=false, slice=50; +makethin, sequence=myseq, style=teapot; + +!save, sequence= myseq,file= "myseq_thin.seq"; +!call, file= "myseq_thin.seq"; + +use, sequence=myseq; +track, onepass, dump, file= "madx-thin"; +start, x=1e-3, y=0, px=0, py=0, t=0, pt=-1e-3; +start, x=1e-3, y=0, px=0, py=0, t=0, pt=0.0; +start, x=1e-3, y=0, px=0, py=0, t=0, pt=1e-3; +run, turns=1; +endtrack; +stop; + diff --git a/tests/test-thick-solenoid/test-thick-solenoid.ref b/tests/test-thick-solenoid/test-thick-solenoid.ref new file mode 100644 index 000000000..7c71aa9db --- /dev/null +++ b/tests/test-thick-solenoid/test-thick-solenoid.ref @@ -0,0 +1,77 @@ + + ++++++++++++++++++++++++++++++++++++++++++++ + + MAD-X 5.03.07 (64 bit, Darwin) + + + Support: mad@cern.ch, http://cern.ch/mad + + + Release date: 2017.10.20 + + + Execution date: 2017.11.06 18:38:10 + + ++++++++++++++++++++++++++++++++++++++++++++ +option, -debug, -echo; + +enter TRACK module +one pass is on + +++++++ table: tracksumm + + number turn x px + 1 0 0.001 0 + 2 0 0.001 0 + 3 0 0.001 0 + 1 1 -0.0004179664169 -9.093718502e-05 + 2 1 -0.0004161468365 -9.092974268e-05 + 3 1 -0.0004143292288 -9.092178557e-05 + + y py t pt + 0 0 0 -0.001 + 0 0 0 2.220446049e-16 + 0 0 0 0.001 + 0 0 -6.495195848e-08 -0.001 + 0 3.388131789e-21 -5.946003867e-08 2.220446049e-16 + 6.783039842e-20 0 -5.398501379e-08 0.001 + + s e + 0 0 + 0 0 + 0 0 + 20 0 + 20 0 + 20 0 +exit TRACK module + +makethin: style chosen : teapot +makethin: slicing sequence : myseq +enter TRACK module +one pass is on + +++++++ table: tracksumm + + number turn x px + 1 0 0.001 0 + 2 0 0.001 0 + 3 0 0.001 0 + 1 1 -0.0004179968082 -9.084288303e-05 + 2 1 -0.0004161771647 -9.092639127e-05 + 3 1 -0.0004143594937 -9.100936927e-05 + + y py t pt + 0 0 0 -0.001 + 0 0 0 2.220446049e-16 + 0 0 0 0.001 + -2.711952481e-19 -3.599890026e-20 -1.054344426e-07 -0.001 + -2.050193423e-19 5.145725155e-20 -1.000038811e-07 2.220446049e-16 + -4.805886845e-19 -2.922263668e-20 -9.458958447e-08 0.001 + + s e + 0 0 + 0 0 + 0 0 + 20 0 + 20 0 + 20 0 +exit TRACK module + + + Number of warnings: 0 + + ++++++++++++++++++++++++++++++++++++++++++++ + + MAD-X finished normally + + ++++++++++++++++++++++++++++++++++++++++++++