Bug Fix: Wiggler diffusion matrix (#759)

* Adding wiggler matrix diffusion code and modifying wiggler integrators to initialize z position for every particle

* Silencing warnings for .mexw64 and .pyd extensions

* Defining energy in mex file from Param instead of ElemData

* Wiggler difussion matrices added to ohmienvelope

* Matlab/Python mex function update

* update matlab & python runners

* Merged gwigrad2 to gwigR

* GwigInit function update

* Adding rotation and translation to integrators

* Spelling typos

* Wiggler index deletion enhanced

* File modification log update

* Updated radiation kicks computation and documentation

* Minor changes to documentation

* Updating rotation and translation diffusion

* Minor changes to integrator file

* Translation diffusion deleted

* radiation kicks back to original formula

---------

Co-authored-by: ZeusMarti <[email protected]>
Co-authored-by: oscarxblanco <[email protected]>

.gitignore

@@ -5,9 +5,11 @@
 *.mex
 *.mexa64
 *.mexmaci64
+*.mexw64
 
 # Compiled python
 *.pyc
+*.pyd
 
 # pytest cache
 .cache

atintegrators/GWigSymplecticPass.c

@@ -6,6 +6,11 @@
  *---------------------------------------------------------------------------
  * Modification Log:
  * -----------------
+ * .03  2024-05-06     J. Arenillas, ALBA, [email protected]
+ *              Adding rotations and translations to wiggler.
+ *				Bug fix in wiggler initialisation.
+ *				Energy parameter bug fix.
+ *				
  * .02  2003-06-18     J. Li
  *				Cleanup the code
  *
@@ -123,11 +128,13 @@ void GWigSymplecticPass(double *r, double Energy, double Ltot, double Lw,
     /* Energy is defined in the lattice in eV but GeV is used by the gwig code. */
     Energy = Energy / 1e9;
 
-    GWigInit(&Wig, Energy, Ltot, Lw, Bmax, Nstep, Nmeth, NHharm, NVharm, By, Bx, T1, T2, R1, R2);
-
     for (c = 0;c<num_particles;c++) {
+		GWigInit(&Wig, Energy, Ltot, Lw, Bmax, Nstep, Nmeth, NHharm, NVharm, By, Bx, T1, T2, R1, R2);
         r6 = r+c*6;
         if (!atIsNaN(r6[0])) {
+			/* Misalignment at entrance */
+			if (T1) ATaddvv(r6,T1);
+            if (R1) ATmultmv(r6,R1);
             switch (Nmeth) {
                 case second:
                     GWigPass_2nd(&Wig, r6);
@@ -139,6 +146,9 @@ void GWigSymplecticPass(double *r, double Energy, double Ltot, double Lw,
                     printf("Invalid wiggler integration method %d.\n", Nmeth);
                     break;
             }
+			/* Misalignment at exit */
+            if (R2) ATmultmv(r6,R2);
+            if (T2) ATaddvv(r6,T2);
         }
     }
 }
@@ -158,7 +168,7 @@ ExportMode struct elem *trackFunction(const atElem *ElemData,struct elem *Elem,
         int Nstep, Nmeth;
         int NHharm, NVharm;
 
-        Energy = atGetDouble(ElemData, "Energy"); check_error();
+        Energy=atGetOptionalDouble(ElemData,"Energy",Param->energy); check_error();
         Ltot = atGetDouble(ElemData, "Length"); check_error();
         Lw = atGetDouble(ElemData, "Lw"); check_error();
         Bmax = atGetDouble(ElemData, "Bmax"); check_error();

atintegrators/GWigSymplecticRadPass.c

@@ -6,6 +6,10 @@
  *---------------------------------------------------------------------------
  * Modification Log:
  * -----------------
+ * .03  2024-05-06     J. Arenillas, ALBA, [email protected]
+ *              Adding rotations and translations to wiggler.
+ *				Bug fix in wiggler initialisation.
+ *				Energy parameter bug fix.
  * .02  2003-06-18     J. Li
  *				Cleanup the code
  *
@@ -123,7 +127,7 @@ void GWigInit(struct gwigR *Wig,double design_energy, double Ltot, double Lw,
 
 #define second 2
 #define fourth 4
-void GWigSymplecticRadPass(double *r,double Energy, double Ltot, double Lw,
+void GWigSymplecticRadPass(double *r, double Energy, double Ltot, double Lw,
             double Bmax, int Nstep, int Nmeth, int NHharm, int NVharm,
             double *By, double *Bx, double *T1, double *T2,
             double *R1, double *R2, int num_particles)
@@ -141,11 +145,13 @@ void GWigSymplecticRadPass(double *r,double Energy, double Ltot, double Lw,
     zEndPointV[0] = 0;
     zEndPointV[1] = Ltot;
 
-    GWigInit(&Wig, Energy, Ltot, Lw, Bmax, Nstep, Nmeth, NHharm, NVharm,0, 0, zEndPointH, zEndPointV, By, Bx, T1, T2, R1, R2);
-
     for(c = 0;c<num_particles;c++) {
+		GWigInit(&Wig, Energy, Ltot, Lw, Bmax, Nstep, Nmeth, NHharm, NVharm,0, 0, zEndPointH, zEndPointV, By, Bx, T1, T2, R1, R2);
         r6 = r+c*6;
         if(!atIsNaN(r6[0])) {
+			/* Misalignment at entrance */
+			if (T1) ATaddvv(r6,T1);
+            if (R1) ATmultmv(r6,R1);
             switch (Nmeth) {
                 case second:
                     GWigPass_2nd(&Wig, r6);
@@ -157,6 +163,9 @@ void GWigSymplecticRadPass(double *r,double Energy, double Ltot, double Lw,
                     printf("Invalid wiggler integration method %d.\n", Nmeth);
                     break;
             }
+			/* Misalignment at exit */
+            if (R2) ATmultmv(r6,R2);
+            if (T2) ATaddvv(r6,T2);
         }
     }
 }
@@ -176,7 +185,7 @@ ExportMode struct elem *trackFunction(const atElem *ElemData,struct elem *Elem,
         int Nstep, Nmeth;
         int NHharm, NVharm;
 
-        Energy = atGetDouble(ElemData, "Energy"); check_error();
+        Energy=atGetOptionalDouble(ElemData,"Energy",Param->energy); check_error();
         Ltot = atGetDouble(ElemData, "Length"); check_error();
         Lw = atGetDouble(ElemData, "Lw"); check_error();
         Bmax = atGetDouble(ElemData, "Bmax"); check_error();