Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Remove the dependency of RFCavityPass on harmonic number #319

Merged
merged 3 commits into from
Oct 20, 2021
Merged

Remove the dependency of RFCavityPass on harmonic number #319

Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
d06ab7e
remove the dependence of RFCavityPass on harmonic number
lfarv Oct 18, 2021
62f1529
Small correction
lfarv Oct 18, 2021
926cee3
simplified RFCavityPass
lfarv Oct 19, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
71 changes: 30 additions & 41 deletions atintegrators/RFCavityPass.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,8 +5,10 @@
* Nicola Carmignani
*/

#include <math.h>
#include "atelem.c"
#include "atlalib.c"
#include "driftkickrad.c"

#define TWOPI 6.28318530717959
#define C0 2.99792458e8
Expand All @@ -27,42 +29,31 @@ void RFCavityPass(double *r_in, double le, double nv, double freq, double h, dou
r is a 6-by-N matrix of initial conditions reshaped into
1-d array of 6*N elements
*/
{ int c, c6;
double halflength , p_norm, NormL;
/* I get T0 and nturn from matlab global variables */
/*T0=getT0FromMatlab();*/

if(le == 0)
{
for(c = 0;c<num_particles;c++)
{ c6 = c*6;
if(!atIsNaN(r_in[c6]))
r_in[c6+4] += -nv*sin(TWOPI*freq*((r_in[c6+5]-lag)/C0 - (h/freq-T0)*nturn ));
}
{
int c;

if (le == 0) {
for (c = 0; c<num_particles; c++) {
double *r6 = r_in+c*6;
if(!atIsNaN(r6[0]))
r6[4] += -nv*sin(TWOPI*freq*((r6[5]-lag)/C0 - (h/freq-T0)*nturn));
}
}
else
{ halflength = le/2;
for(c = 0;c<num_particles;c++)
{ c6 = c*6;
if(!atIsNaN(r_in[c6]))
{ p_norm = 1/(1+r_in[c6+4]);
NormL = halflength*p_norm;
/* Propagate through a drift equal to half cavity length */
r_in[c6+0]+= NormL*r_in[c6+1];
r_in[c6+2]+= NormL*r_in[c6+3];
r_in[c6+5]+= NormL*p_norm*(r_in[c6+1]*r_in[c6+1]+r_in[c6+3]*r_in[c6+3])/2;
/* Longitudinal momentum kick */
r_in[c6+4] += -nv*sin(TWOPI*freq*((r_in[c6+5]-lag)/C0 - (h/freq-T0)*nturn ));
p_norm = 1/(1+r_in[c6+4]);
NormL = halflength*p_norm;
/* Propagate through a drift equal to half cavity length */
r_in[c6+0]+= NormL*r_in[c6+1];
r_in[c6+2]+= NormL*r_in[c6+3];
r_in[c6+5]+= NormL*p_norm*(r_in[c6+1]*r_in[c6+1]+r_in[c6+3]*r_in[c6+3])/2;
}
}
else {
double halflength = le/2;
for (c = 0;c<num_particles;c++) {
double *r6 = r_in+c*6;
if(!atIsNaN(r6[0])) {
/* Propagate through a drift equal to half cavity length */
drift6(r6, halflength);
/* Longitudinal momentum kick */
r6[4] += -nv*sin(TWOPI*freq*((r6[5]-lag)/C0 - (h/freq-T0)*nturn));
/* Propagate through a drift equal to half cavity length */
drift6(r6, halflength);
}
}
}
}
}

#if defined(MATLAB_MEX_FILE) || defined(PYAT)
ExportMode struct elem *trackFunction(const atElem *ElemData,struct elem *Elem,
Expand All @@ -71,19 +62,18 @@ ExportMode struct elem *trackFunction(const atElem *ElemData,struct elem *Elem,
int nturn=Param->nturn;
double T0=Param->T0;
if (!Elem) {
double Length, Voltage, Energy, Frequency, HarmNumber, TimeLag;
double Length, Voltage, Energy, Frequency, TimeLag;
Length=atGetDouble(ElemData,"Length"); check_error();
Voltage=atGetDouble(ElemData,"Voltage"); check_error();
Energy=atGetDouble(ElemData,"Energy"); check_error();
Frequency=atGetDouble(ElemData,"Frequency"); check_error();
HarmNumber=atGetDouble(ElemData,"HarmNumber"); check_error();
TimeLag=atGetOptionalDouble(ElemData,"TimeLag",0); check_error();
Elem = (struct elem*)atMalloc(sizeof(struct elem));
Elem->Length=Length;
Elem->Voltage=Voltage;
Elem->Energy=Energy;
Elem->Frequency=Frequency;
Elem->HarmNumber=HarmNumber;
Elem->HarmNumber=round(Frequency*T0);
Elem->TimeLag=TimeLag;
}
RFCavityPass(r_in, Elem->Length, Elem->Voltage/Elem->Energy, Elem->Frequency, Elem->HarmNumber, Elem->TimeLag, nturn, T0, num_particles);
Expand All @@ -107,9 +97,9 @@ void mexFunction( int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
double Voltage=atGetDouble(ElemData,"Voltage");
double Energy=atGetDouble(ElemData,"Energy");
double Frequency=atGetDouble(ElemData,"Frequency");
double HarmNumber=atGetDouble(ElemData,"HarmNumber");
double TimeLag=atGetOptionalDouble(ElemData,"TimeLag",0);
double T0=HarmNumber/Frequency;
double T0=1.0/Frequency; /* Does not matter since nturns == 0 */
double HarmNumber=round(Frequency*T0)
if (mxGetM(prhs[1]) != 6) mexErrMsgIdAndTxt("AT:WrongArg","Second argument must be a 6 x N matrix");
/* ALLOCATE memory for the output array of the same size as the input */
plhs[0] = mxDuplicateArray(prhs[1]);
Expand All @@ -119,12 +109,11 @@ void mexFunction( int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
}
else if (nrhs == 0)
{ /* return list of required fields */
plhs[0] = mxCreateCellMatrix(5,1);
plhs[0] = mxCreateCellMatrix(4,1);
mxSetCell(plhs[0],0,mxCreateString("Length"));
mxSetCell(plhs[0],1,mxCreateString("Voltage"));
mxSetCell(plhs[0],2,mxCreateString("Energy"));
mxSetCell(plhs[0],3,mxCreateString("Frequency"));
mxSetCell(plhs[0],4,mxCreateString("HarmNumber"));
if(nlhs>1) /* optional fields */
{
plhs[1] = mxCreateCellMatrix(1,1);
Expand Down