Option, ECHO=FALSE;
Option, PSDUMPFREQ=100000;
Option, SPTDUMPFREQ = 10;
Option, PSDUMPEACHTURN=false;
Option, PSDUMPLOCALFRAME=FALSE;
Option, REPARTFREQ=20;
Option, ECHO=FALSE;
Option, STATDUMPFREQ=1;
Option, CZERO=FALSE;
Option, TELL=TRUE;
Option, VERSION=10900;

Title,string="OPAL-cycl: the first turn acceleration in PSI 590MeV Ring";

REAL Edes=.072;
REAL gamma=(Edes+PMASS)/PMASS;
REAL beta=sqrt(1-(1/gamma^2));
REAL gambet=gamma*beta;
REAL P0 = gamma*beta*PMASS;
REAL brho = (PMASS*1.0e9*gambet) / CLIGHT;

//value,{gamma,brho,Edes,beta,gambet};

REAL phi01=139.4281;
REAL phi02=phi01+180.0;
REAL phi04=phi01;
REAL phi05=phi01+180.0;
REAL phi03=phi01+10.0;

REAL volt1st=0.9;
REAL volt3rd=0.9*4.0*0.112;

REAL turns = 1;
REAL nstep=2000;

REAL frequency=50.650;
REAL frequency3=3.0*frequency;

ring: Cyclotron, TYPE="RING", CYHARMON=6, PHIINIT=0.0,
	PRINIT=-0.000174, RINIT=2130.0, SYMMETRY=8.0, RFFREQ=frequency,
	FMAPFN="s03av.nar",
	TCR1V={2022.09, 2012.09, 2002.09, 2293.59, 2501.36, 2694.86, 2893.45, 3092.15,
	       3293.29, 3491.55, 3690.11, 3892.46, 4060.34, 4217.20, 4350, 4471},
	TCR2V={2132.09, 2212.09, 2292.09, 2498.09, 2694.36, 2893.36, 3091.95, 3290.65,
	       3491.29, 3689.55, 3888.11, 4058.46, 4213.34, 4349.20, 4470, 4600},
	MBTCV={0.0e-3, 0.0e-3, 0.0e-3, 0.0e-3, 2.0e-3, 0.0e-3, 0.0e-3, 0.0e-3,
	       0.0e-3, 0.0e-3, 0.0e-3, 0.0e-3, 0.0e-3, 0.0e-3, 0.0e-3, 0.0e-3},
	SLPTCV={1, 1, 1, 6, 6, 6, 6, 6,
	        6, 6, 6, 6, 6, 6, 6, 6};

rf0: RFCavity, VOLT=volt1st, FMAPFN="Cav1.dat", TYPE="SINGLEGAP",
	FREQ=frequency, RMIN = 1900.0, RMAX = 4500.0, ANGLE=35.0,  PDIS = 416.0,
	GAPWIDTH = 220.0, PHI0=phi01;

rf1: RFCavity, VOLT=volt1st, FMAPFN="Cav1.dat", TYPE="SINGLEGAP",
	FREQ=frequency, RMIN = 1900.0, RMAX = 4500.0, ANGLE=125.0, PDIS = 416.0,
	GAPWIDTH = 220.0, PHI0=phi02;

rf2: RFCavity, VOLT=volt3rd, FMAPFN="Cav3.dat", TYPE="SINGLEGAP",
	FREQ=frequency3,RMIN = 1900.0, RMAX = 4500.0, ANGLE=170.0, PDIS = 452.0,
	GAPWIDTH = 250.0, PHI0=phi03;

rf3: RFCavity, VOLT=volt1st, FMAPFN="Cav1.dat", TYPE="SINGLEGAP",
	FREQ=frequency, RMIN = 1900.0, RMAX = 4500.0, ANGLE=215.0, PDIS = 416.0,
	GAPWIDTH = 220.0, PHI0=phi04;

rf4: RFCavity, VOLT=volt1st, FMAPFN="Cav1.dat", TYPE="SINGLEGAP",
	FREQ=frequency, RMIN = 1900.0, RMAX = 4500.0, ANGLE=305.0, PDIS = 416.0,
	GAPWIDTH = 220.0, PHI0=phi05;

l1:   Line = (ring,rf0,rf1,rf2,rf3,rf4);

Dist1:DISTRIBUTION, TYPE=gauss,
	sigmax = 2.0e-03,
	sigmapx = 1.0e-7,
	corrx = 0.0,
	sigmay = 2.0e-03,
	sigmapy = 1.0e-7,
	corry = 0.0,
	sigmat = 2.0e-03,
	sigmapt = 3.394e-4,
	corrt=0.0;

Fs1:FIELDSOLVER, FSTYPE=FFT, MX=32, MY=32, MT=16,
	PARFFTX=true, PARFFTY=true, PARFFTT=true,
	BCFFTX=open, BCFFTY=open, BCFFTT=open, BBOXINCR=2;

beam1: BEAM, PARTICLE=PROTON, pc=P0, SPACECHARGE=true, NPART=1e5, BCURRENT=1.0E-3, CHARGE=1.0, BFREQ= frequency;

Select, Line=l1;

TRACK,LINE=l1, BEAM=beam1, MAXSTEPS=nstep*turns, STEPSPERTURN=nstep,TIMEINTEGRATOR="RK-4";
 
 run, method = "CYCLOTRON-T", beam=beam1, fieldsolver=Fs1, distribution=Dist1;

endtrack;

Stop;
