Newer
Older
#define ORDER 1
int no_tps = ORDER; // arbitrary TPSA order is defined locally
extern bool freq_map;
#include "tracy_lib.h"
//***************************************************************************************
//
// MAIN CODE
//
//****************************************************************************************
int main(int argc, char *argv[])
{
const long seed = 1121;
iniranf(seed); setrancut(2.0);
// turn on globval.Cavity_on and globval.radiation to get proper synchr radiation damping
// IDs accounted too if: wiggler model and symplectic integrator (method = 1)
globval.H_exact = false;
globval.quad_fringe = false; // quadrupole fringe field
globval.radiation = false; // synchrotron radiation
globval.emittance = false; // emittance
// const double x_max_FMA = 10e-3, delta_FMA = 10e-2;
// const int n_x = 801, n_dp = 80, n_tr = 2048;
bool chroma;
double dP = 0.0;
long lastpos = -1L;
char str1[S_SIZE];
// for time handling
uint32_t start, stop;
/************************************************************************
start read in files and flags
*************************************************************************/
if (argc > 1){
read_script(argv[1], true);}
else{
fprintf(stdout, "Not enough parameters\nSyntax is program parameterfile\n");
return 1;
}
/************************************************************************
end read in files and flags
*************************************************************************/
prtmfile("flat_file.dat"); // writes flat file /* very important file for debug*/
printlatt(); /* SOLEIL print out lattice functions */
printglob();
// Flag factory
// Chamber factory
if (ChamberFlag == false)
ChamberOff(); // turn off vacuum chamber setting, use the default one
else if (ChamberNoU20Flag == true)
DefineChNoU20(); // using vacuum chamber setting but without undulator U20
else if (ReadChamberFlag == true)
ReadCh(chamber_file); /* read vacuum chamber from a file "Apertures.dat" , soleil version*/
PrintCh(); // print chamber into chamber.out
//get_matching_params_scl(); // get tunes and beta functions at entrance
get_alphac2(); // compute up to 3rd order mcf
//cout << endl << "computing tune shifts" << endl;
//dnu_dA(10e-3, 5e-3, 0.002);
//get_ksi2(0.0); // this gets the chromas and writes them into chrom2.out
// compute tunes by tracking (should be the same as by DA)
if (TuneTracFlag == true) {
GetTuneTrac(1026L, 0.0, &nux, &nuz);
fprintf(stdout,"From tracking: nux = % f nuz = % f \n",nux,nuz);
}
// compute chromaticities by tracking (should be the same as by DA)
if (ChromTracFlag == true){
start = stampstart();
GetChromTrac(2L, 1026L, 1e-5, &ksix, &ksiz);
stop = stampstop(start);
fprintf(stdout,"From tracking: ksix= % f ksiz= % f \n",ksix,ksiz);
}
if (FitTuneFlag == true){
fprintf(stdout, "\n Fitting tunes\n");
FitTune(ElemIndex("qp7"),ElemIndex("qp9"), targetnux, targetnuz);
Ring_GetTwiss(chroma=true, 0.0); /* Compute and get Twiss parameters */
printglob(); /* print parameter list */
}
if (FitChromFlag == true){
fprintf(stdout, "\n Fitting chromaticities\n");
FitChrom(ElemIndex("sx9"),ElemIndex("sx10"), targetksix, targetksiz);
Ring_GetTwiss(chroma=true, 0.0); /* Compute and get Twiss parameters */
printglob(); /* print parameter list */
}
// coupling calculation
if (CouplingFlag == true){
Ring_GetTwiss(chroma=true, 0.0); /* Compute and get Twiss parameters */
printlatt(); /* dump linear lattice functions into "linlat.dat" */
GetEmittance(ElemIndex("cav"), true);
Coupling_Edwards_Teng();
printglob(); /* print parameter list */
}
// add coupling by random rotating of the quadrupoles
if (ErrorCouplingFlag == true){
SetErr();
Ring_GetTwiss(chroma=true, 0.0); /* Compute and get Twiss parameters */
printlatt(); /* dump linear lattice functions into "linlat.dat" */
Coupling_Edwards_Teng();
GetEmittance(ElemIndex("cav"), true);
printglob(); /* print parameter list */
}
// WARNING Fit tunes and chromaticities before applying errors !!!!
//set multipoles in all magnets
if (MultipoleFlag == true ){
if (ThinsextFlag ==true){
fprintf(stdout, "\n Setting Multipoles for lattice with thin sextupoles \n");
Multipole_thinsext(fic_hcorr,fic_vcorr,fic_skew); /* for thin sextupoles */
Ring_GetTwiss(chroma=true, 0.0); /* Compute and get Twiss parameters */
printglob();
}
else{
fprintf(stdout, "\n Setting Multipoles for lattice with thick sextupoles \n");
Multipole_thicksext(fic_hcorr,fic_vcorr,fic_skew); /* for thick sextupoles */
Ring_GetTwiss(chroma=true, 0.0); /* Compute and get Twiss parameters */
printglob();
}
}
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
/******************************************************************************************/
// COMPUTATION PART after setting the model
/******************************************************************************************/
//first print the full lattice with error as a flat file
prtmfile("flat_file_error.dat"); // writes flat file /* very important file for debug*/
// computes TuneShift with amplitudes
if (AmplitudeTuneShiftFlag == true){
if (ChamberFlag == true ){
NuDx(_AmplitudeTuneShift_nxpoint,
_AmplitudeTuneShift_nypoint, _AmplitudeTuneShift_nturn,
_AmplitudeTuneShift_xmax, _AmplitudeTuneShift_ymax,
_AmplitudeTuneShift_delta);
//NuDx(31L,21L,516L,0.025,0.005,dP);
}
else{ // Utility ?
NuDx(50L,30L,516L,0.035,0.02,dP);
}
}
if (EnergyTuneShiftFlag == true){
if (ChamberFlag == true ){
NuDp(_EnergyTuneShift_npoint,
_EnergyTuneShift_nturn, _EnergyTuneShift_deltamax);
//NuDp(31L,1026L,0.06);
}
else{ // utility ?
NuDp(31L,1026L,0.06);
}
}
// Computes FMA
if (FmapFlag == true){
if (ChamberFlag == true ){
if (ExperimentFMAFlag == true)
fmap( _FmapFlag_nxpoint, _FmapFlag_nypoint,
_FmapFlag_nturn, _FmapFlag_xmax, _FmapFlag_ymax,
_FmapFlag_delta, _FmapFlag_diffusion);
//fmap(40,12,258,-20e-3,5e-3,0.0,true); // for experimental
fmap( _FmapFlag_nxpoint, _FmapFlag_nypoint,
_FmapFlag_nturn, _FmapFlag_xmax, _FmapFlag_ymax,
_FmapFlag_delta, _FmapFlag_diffusion);
// fmap(100,50,1026,20e-3,5e-3,0.0,true);
fmap( _FmapFlag_nxpoint, _FmapFlag_nypoint,
_FmapFlag_nturn, _FmapFlag_xmax, _FmapFlag_ymax,
_FmapFlag_delta, _FmapFlag_diffusion);
fmap( _FmapFlag_nxpoint, _FmapFlag_nypoint,
_FmapFlag_nturn, _FmapFlag_xmax, _FmapFlag_ymax,
_FmapFlag_delta, _FmapFlag_diffusion);
// fmap(200,100,1026,32e-3,7e-3,0.0,true);
}
}
if (CodeComparaisonFlag){
fmap(200,100,1026,-32e-3,7e-3,0.0,true);
}
MomentumAcceptance(
_MomentumAccFlag_istart, _MomentumAccFlag_istop,
_MomentumAccFlag_deltaminp, _MomentumAccFlag_deltamaxp, _MomentumAccFlag_nstepp,
_MomentumAccFlag_deltaminn, _MomentumAccFlag_deltamaxn, _MomentumAccFlag_nstepn);
// MomentumAcceptance(1L, 108L, 0.01, 0.05, 100L, -0.01, -0.05, 100L);
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
}
// induced amplitude
if (InducedAmplitudeFlag == true){
InducedAmplitude(193L);
}
if (EtaFlag == true){
// compute cod and twiss parameters for different energy offsets
for (int ii=0; ii<=40; ii++) {
dP = -0.02+ 0.001*ii;
Ring_GetTwiss(chroma=false, dP); /* Compute and get Twiss parameters */
printlatt(); /* dump linear lattice functions into "linlat.dat" */
getcod (dP, lastpos);
// printcod();
prt_cod("cod.out", globval.bpm, true);
//system("mv linlat.out linlat_ooo.out");
sprintf(str1, "mv cod.out cod_%02d.out", ii);
system(str1);
sprintf(str1, "mv linlat.out linlat_%02d.out", ii);
system(str1);
}
}
if (PhaseSpaceFlag == true){
start = stampstart();
Phase(0.001,0,0.001,0,0.001,0, 1);
printf("the simulation time for phase space in tracy 3 is \n");
stop = stampstop(start);