fix large beam tilt scaling problem

Author: Ivanlh20

src/cgpu_fcns.cuh

@@ -1867,12 +1867,13 @@ namespace mt
 		template <class TGrid, class TVector_c>
 		DEVICE_CALLABLE FORCE_INLINE 
 		void exp_r_factor_2d(const int &ix, const int &iy, const TGrid &grid_2d, 
-		const Value_type<TGrid> &gx, const Value_type<TGrid> &gy, TVector_c &psi_i, TVector_c &psi_o)
+		const Value_type<TGrid> &gx, const Value_type<TGrid> &gy, TVector_c &psi_i, TVector_c &psi_o, bool scaling)
 		{
 			const int ixy = grid_2d.ind_col(ix, iy);
 			const auto Rx = grid_2d.Rx_shift(ix)-grid_2d.Rx_c();
 			const auto Ry = grid_2d.Ry_shift(iy)-grid_2d.Ry_c();
-			psi_o[ixy] = psi_i[ixy]*euler(gx*Rx + gy*Ry)/grid_2d.nxy_r();
+			const auto scale = (scaling)?grid_2d.nxy_r():1;
+			psi_o[ixy] = psi_i[ixy]*euler(gx*Rx + gy*Ry)/scale;
 		}
 
 		template <class TGrid, class TVector_r, class TVector_c>

src/cpu_fcns.hpp

@@ -921,11 +921,11 @@ namespace mt
 	template <class TGrid, class TVector_c>
 	enable_if_host_vector<TVector_c, void>
 		exp_r_factor_2d(Stream<e_host> &stream, TGrid &grid_2d, Value_type<TGrid> gx, Value_type<TGrid> gy,
-			TVector_c &fPsi_i, TVector_c &fPsi_o)
+			TVector_c &fPsi_i, TVector_c &fPsi_o, bool scaling=true)
 	{
 		stream.set_n_act_stream(grid_2d.nx);
 		stream.set_grid(grid_2d.nx, grid_2d.ny);
-		stream.exec_matrix(host_device_detail::exp_r_factor_2d<TGrid, TVector_c>, grid_2d, gx, gy, fPsi_i, fPsi_o);
+		stream.exec_matrix(host_device_detail::exp_r_factor_2d<TGrid, TVector_c>, grid_2d, gx, gy, fPsi_i, fPsi_o, scaling);
 	}
 
 	template <class TGrid, class TVector_c>

src/gpu_fcns.cuh

@@ -649,14 +649,14 @@ namespace mt
 		// phase factor 2d
 		template <class TGrid, class T>
 		__global__ void exp_r_factor_2d(TGrid grid_2d, Value_type<TGrid> gx, 
-		Value_type<TGrid> gy, rVector<T> psi_i, rVector<T> psi_o)
+		Value_type<TGrid> gy, rVector<T> psi_i, rVector<T> psi_o, scaling)
 		{
 			int iy = threadIdx.x + blockIdx.x*blockDim.x;
 			int ix = threadIdx.y + blockIdx.y*blockDim.y;
 
 			if((ix < grid_2d.nx) && (iy < grid_2d.ny))
 			{
-				host_device_detail::exp_r_factor_2d(ix, iy, grid_2d, gx, gy, psi_i, psi_o);
+				host_device_detail::exp_r_factor_2d(ix, iy, grid_2d, gx, gy, psi_i, psi_o, scaling);
 			}
 		}
 
@@ -1325,11 +1325,11 @@ namespace mt
 
 	template <class TGrid, class TVector_c>
 	enable_if_device_vector<TVector_c, void>
-	exp_r_factor_2d(Stream<e_device> &stream, TGrid &grid_2d, Value_type<TGrid> gx, Value_type<TGrid> gy, TVector_c &fPsi_i, TVector_c &fPsi_o)
+	exp_r_factor_2d(Stream<e_device> &stream, TGrid &grid_2d, Value_type<TGrid> gx, Value_type<TGrid> gy, TVector_c &fPsi_i, TVector_c &fPsi_o, bool scaling=true)
 	{
 		auto grid_bt = grid_2d.cuda_grid();
 
-		device_detail::exp_r_factor_2d<TGrid, typename TVector_c::value_type><<<grid_bt.Blk, grid_bt.Thr>>>(grid_2d, gx, gy, fPsi_i, fPsi_o);
+		device_detail::exp_r_factor_2d<TGrid, typename TVector_c::value_type><<<grid_bt.Blk, grid_bt.Thr>>>(grid_2d, gx, gy, fPsi_i, fPsi_o, scaling);
 	}
 
 	template <class TGrid, class TVector_c>

src/wave_function.cuh

@@ -74,7 +74,7 @@ namespace mt
 				Transmission_Function<T, dev>::set_input_data(input_multislice_i, stream_i, fft2_i);
 			}
 
-			void phase_multiplication(const T_r &gxu, const T_r &gyu, TVector_c &psi_i, TVector_c &psi_o)
+			void phase_multiplication(const T_r &gxu, const T_r &gyu, TVector_c &psi_i, TVector_c &psi_o, bool scaling=true)
 			{
 				if(this->input_multislice->dp_Shift || isZero(gxu, gyu))
 				{
@@ -85,19 +85,20 @@ namespace mt
 					return;
 				}
 
-				mt::exp_r_factor_2d(*(this->stream), this->input_multislice->grid_2d, c_2Pi*gxu, c_2Pi*gyu, psi_i, psi_o);
+				mt::exp_r_factor_2d(*(this->stream), this->input_multislice->grid_2d, c_2Pi*gxu, c_2Pi*gyu, psi_i, psi_o, scaling);
 			}
 
-			void phase_multiplication(const T_r &gxu, const T_r &gyu, TVector_c &psi_io)
+			void phase_multiplication(const T_r &gxu, const T_r &gyu, TVector_c &psi_io, bool scaling=true)
 			{
-				phase_multiplication(gxu, gyu, psi_io, psi_io);
+				phase_multiplication(gxu, gyu, psi_io, psi_io, scaling);
 			}
 
 			TVector_c* get_psi(const eSpace &space, const T_r &gxu, const T_r &gyu, 
 			T_r z, TVector_c &psi_i)
 			{
 				TVector_c *psi_o = &(this->trans_0);
-				phase_multiplication(gxu, gyu, psi_i, *psi_o);
+				// real space not need to include scaling and phase multiplication oposite sign as the propagation
+				phase_multiplication(-gxu, -gyu, psi_i, *psi_o, false); 
 				propagator(space, gxu, gyu, z, *psi_o);
 
 				return psi_o;