-
Notifications
You must be signed in to change notification settings - Fork 20
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.
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
Element to apply a simple x-y rotation #747
Open
cemitch99
wants to merge
5
commits into
ECP-WarpX:development
Choose a base branch
from
cemitch99:add_plane_rotation
base: development
Could not load branches
Branch not found: {{ refName }}
Loading
Could not load tags
Nothing to show
Loading
Are you sure you want to change the base?
Some commits from the old base branch may be removed from the timeline,
and old review comments may become outdated.
Open
Changes from 2 commits
Commits
Show all changes
5 commits
Select commit
Hold shift + click to select a range
5d13c90
Add element for simple x-y rotation.
cemitch99 ac8fd76
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] b0a8c51
Add documentation.
cemitch99 211ad4e
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] 44d1b9e
Add ParmParse for element intialization.
cemitch99 File filter
Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,98 @@ | ||
#!/usr/bin/env python3 | ||
# | ||
# Copyright 2022-2023 ImpactX contributors | ||
# Authors: Axel Huebl, Chad Mitchell | ||
# License: BSD-3-Clause-LBNL | ||
# | ||
|
||
|
||
import numpy as np | ||
import openpmd_api as io | ||
from scipy.stats import moment | ||
|
||
|
||
def get_moments(beam): | ||
"""Calculate standard deviations of beam position & momenta | ||
and emittance values | ||
|
||
Returns | ||
------- | ||
sigx, sigy, sigt, emittance_x, emittance_y, emittance_t | ||
""" | ||
sigx = moment(beam["position_x"], moment=2) ** 0.5 # variance -> std dev. | ||
sigpx = moment(beam["momentum_x"], moment=2) ** 0.5 | ||
sigy = moment(beam["position_y"], moment=2) ** 0.5 | ||
sigpy = moment(beam["momentum_y"], moment=2) ** 0.5 | ||
sigt = moment(beam["position_t"], moment=2) ** 0.5 | ||
sigpt = moment(beam["momentum_t"], moment=2) ** 0.5 | ||
|
||
epstrms = beam.cov(ddof=0) | ||
emittance_x = (sigx**2 * sigpx**2 - epstrms["position_x"]["momentum_x"] ** 2) ** 0.5 | ||
emittance_y = (sigy**2 * sigpy**2 - epstrms["position_y"]["momentum_y"] ** 2) ** 0.5 | ||
emittance_t = (sigt**2 * sigpt**2 - epstrms["position_t"]["momentum_t"] ** 2) ** 0.5 | ||
|
||
return (sigx, sigy, sigt, emittance_x, emittance_y, emittance_t) | ||
|
||
|
||
# initial/final beam | ||
series = io.Series("diags/openPMD/monitor.h5", io.Access.read_only) | ||
last_step = list(series.iterations)[-1] | ||
initial = series.iterations[1].particles["beam"].to_df() | ||
final = series.iterations[last_step].particles["beam"].to_df() | ||
|
||
# compare number of particles | ||
num_particles = 10000 | ||
assert num_particles == len(initial) | ||
assert num_particles == len(final) | ||
|
||
print("Initial Beam:") | ||
sigx, sigy, sigt, emittance_x, emittance_y, emittance_t = get_moments(initial) | ||
print(f" sigx={sigx:e} sigy={sigy:e} sigt={sigt:e}") | ||
print( | ||
f" emittance_x={emittance_x:e} emittance_y={emittance_y:e} emittance_t={emittance_t:e}" | ||
) | ||
|
||
atol = 0.0 # ignored | ||
rtol = 1.5e12 * num_particles**-0.5 # from random sampling of a smooth distribution | ||
print(f" rtol={rtol} (ignored: atol~={atol})") | ||
|
||
assert np.allclose( | ||
[sigx, sigy, sigt, emittance_x, emittance_y, emittance_t], | ||
[ | ||
4.0e-03, | ||
2.0e-03, | ||
1.0e-03, | ||
8.0e-07, | ||
1.0e-06, | ||
2.0e-06, | ||
], | ||
rtol=rtol, | ||
atol=atol, | ||
) | ||
|
||
|
||
print("") | ||
print("Final Beam:") | ||
sigx, sigy, sigt, emittance_x, emittance_y, emittance_t = get_moments(final) | ||
print(f" sigx={sigx:e} sigy={sigy:e} sigt={sigt:e}") | ||
print( | ||
f" emittance_x={emittance_x:e} emittance_y={emittance_y:e} emittance_t={emittance_t:e}" | ||
) | ||
|
||
atol = 0.0 # ignored | ||
rtol = 1.5 * num_particles**-0.5 # from random sampling of a smooth distribution | ||
print(f" rtol={rtol} (ignored: atol~={atol})") | ||
|
||
assert np.allclose( | ||
[sigx, sigy, sigt, emittance_x, emittance_y, emittance_t], | ||
[ | ||
2.0e-03, | ||
4.0e-03, | ||
1.0e-03, | ||
1.0e-06, | ||
8.0e-07, | ||
2.0e-06, | ||
], | ||
rtol=rtol, | ||
atol=atol, | ||
) |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,33 @@ | ||
############################################################################### | ||
# Particle Beam(s) | ||
############################################################################### | ||
beam.npart = 10000 | ||
beam.units = static | ||
beam.kin_energy = 2.0e3 | ||
beam.charge = 1.0e-9 | ||
beam.particle = electron | ||
beam.distribution = waterbag | ||
beam.lambdaX = 4.0e-3 | ||
beam.lambdaY = 2.0e-3 | ||
beam.lambdaT = 1.0e-3 | ||
beam.lambdaPx = 2.0e-4 | ||
beam.lambdaPy = 5.0e-4 | ||
beam.lambdaPt = 2.0e-3 | ||
|
||
|
||
############################################################################### | ||
# Beamline: lattice elements and segments | ||
############################################################################### | ||
lattice.elements = monitor rotation1 monitor | ||
|
||
monitor.type = beam_monitor | ||
monitor.backend = h5 | ||
|
||
rotation1.type = plane_xyrotation | ||
rotation1.angle = 90.0 | ||
|
||
############################################################################### | ||
# Algorithms | ||
############################################################################### | ||
algo.particle_shape = 2 | ||
algo.space_charge = false |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,58 @@ | ||
#!/usr/bin/env python3 | ||
# | ||
# Copyright 2022-2023 ImpactX contributors | ||
# Authors: Ryan Sandberg, Axel Huebl, Chad Mitchell | ||
# License: BSD-3-Clause-LBNL | ||
# | ||
# -*- coding: utf-8 -*- | ||
|
||
from impactx import ImpactX, distribution, elements | ||
|
||
sim = ImpactX() | ||
|
||
# set numerical parameters and IO control | ||
sim.particle_shape = 2 # B-spline order | ||
sim.space_charge = False | ||
# sim.diagnostics = False # benchmarking | ||
sim.slice_step_diagnostics = True | ||
|
||
# domain decomposition & space charge mesh | ||
sim.init_grids() | ||
|
||
# load a 2 GeV electron beam with an initial | ||
# unnormalized rms emittance of nm | ||
kin_energy_MeV = 2.0e3 # reference energy | ||
bunch_charge_C = 1.0e-9 # used without space charge | ||
npart = 10000 # number of macro particles | ||
|
||
# reference particle | ||
ref = sim.particle_container().ref_particle() | ||
ref.set_charge_qe(-1.0).set_mass_MeV(0.510998950).set_kin_energy_MeV(kin_energy_MeV) | ||
|
||
# particle bunch | ||
distr = distribution.Waterbag( | ||
lambdaX=4.0e-3, | ||
lambdaY=2.0e-3, | ||
lambdaT=1.0e-3, | ||
lambdaPx=2.0e-4, | ||
lambdaPy=5.0e-4, | ||
lambdaPt=2.0e-3, | ||
) | ||
sim.add_particles(bunch_charge_C, distr, npart) | ||
|
||
# add beam diagnostics | ||
monitor = elements.BeamMonitor("monitor", backend="h5") | ||
|
||
# 90 degree rotation | ||
rotation = elements.PlaneXYRot(name="rotation1", angle=90.0) | ||
|
||
# assign a lattice segment | ||
sim.lattice.append(monitor) | ||
sim.lattice.append(rotation) | ||
sim.lattice.append(monitor) | ||
|
||
# run simulation | ||
sim.track_particles() | ||
|
||
# clean shutdown | ||
sim.finalize() |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,139 @@ | ||
/* Copyright 2022-2023 The Regents of the University of California, through Lawrence | ||
* Berkeley National Laboratory (subject to receipt of any required | ||
* approvals from the U.S. Dept. of Energy). All rights reserved. | ||
* | ||
* This file is part of ImpactX. | ||
* | ||
* Authors: Chad Mitchell, Axel Huebl | ||
* License: BSD-3-Clause-LBNL | ||
*/ | ||
#ifndef IMPACTX_PLANE_XYROT_H | ||
#define IMPACTX_PLANE_XYROT_H | ||
|
||
#include "particles/ImpactXParticleContainer.H" | ||
#include "mixin/alignment.H" | ||
#include "mixin/beamoptic.H" | ||
#include "mixin/thin.H" | ||
#include "mixin/named.H" | ||
#include "mixin/nofinalize.H" | ||
|
||
#include <ablastr/constant.H> | ||
|
||
#include <AMReX_Extension.H> | ||
#include <AMReX_Math.H> | ||
#include <AMReX_REAL.H> | ||
|
||
#include <cmath> | ||
|
||
|
||
namespace impactx | ||
{ | ||
struct PlaneXYRot | ||
: public elements::Named, | ||
public elements::BeamOptic<PlaneXYRot>, | ||
public elements::Thin, | ||
public elements::Alignment, | ||
public elements::NoFinalize | ||
{ | ||
static constexpr auto type = "PlaneXYRot"; | ||
using PType = ImpactXParticleContainer::ParticleType; | ||
|
||
static constexpr amrex::ParticleReal degree2rad = ablastr::constant::math::pi / 180.0; | ||
|
||
/** A simple rotation by an angle phi in the plane transverse to the velocity | ||
* of the reference particle (x-y plane). This is the linear symplectic map | ||
* generated by the longitudinal angular momentum J_z. | ||
* | ||
* @param phi angle of rotation in the x-y plane (about the direction of the ref traj) (degrees) | ||
* @param dx horizontal translation error in m | ||
* @param dy vertical translation error in m | ||
* @param rotation_degree rotation error in the transverse plane [degrees] | ||
* @param name a user defined and not necessarily unique name of the element | ||
*/ | ||
PlaneXYRot ( | ||
amrex::ParticleReal phi, | ||
amrex::ParticleReal dx = 0, | ||
amrex::ParticleReal dy = 0, | ||
amrex::ParticleReal rotation_degree = 0, | ||
std::optional<std::string> name = std::nullopt | ||
) | ||
: Named(name), | ||
Alignment(dx, dy, rotation_degree), | ||
m_phi(phi * degree2rad) | ||
{ | ||
} | ||
|
||
/** Push all particles */ | ||
using BeamOptic::operator(); | ||
|
||
/** This is a prot functor, so that a variable of this type can be used like a | ||
* prot function. | ||
* | ||
* @param x particle position in x | ||
* @param y particle position in y | ||
* @param t particle position in t | ||
* @param px particle momentum in x | ||
* @param py particle momentum in y | ||
* @param pt particle momentum in t | ||
* @param idcpu particle global index (unused) | ||
* @param refpart reference particle | ||
*/ | ||
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE | ||
void operator() ( | ||
amrex::ParticleReal & AMREX_RESTRICT x, | ||
amrex::ParticleReal & AMREX_RESTRICT y, | ||
amrex::ParticleReal & AMREX_RESTRICT t, | ||
amrex::ParticleReal & AMREX_RESTRICT px, | ||
amrex::ParticleReal & AMREX_RESTRICT py, | ||
amrex::ParticleReal & AMREX_RESTRICT pt, | ||
[[maybe_unused]] uint64_t & AMREX_RESTRICT idcpu, | ||
[[maybe_unused]] RefPart const & refpart | ||
) const | ||
{ | ||
using namespace amrex::literals; // for _rt and _prt | ||
|
||
// shift due to alignment errors of the element | ||
shift_in(x, y, px, py); | ||
|
||
// intialize output values | ||
amrex::ParticleReal xout = x; | ||
amrex::ParticleReal yout = y; | ||
amrex::ParticleReal tout = t; | ||
amrex::ParticleReal pxout = px; | ||
amrex::ParticleReal pyout = py; | ||
amrex::ParticleReal ptout = pt; | ||
|
||
// store sin(phi) and cos(phi) | ||
auto const [sin_phi, cos_phi] = amrex::Math::sincos(m_phi); | ||
|
||
// advance position and momentum | ||
xout = x*cos_phi - y*sin_phi; | ||
pxout = px*cos_phi - py*sin_phi; | ||
|
||
yout = x*sin_phi + y*cos_phi; | ||
pyout = px*sin_phi + py*cos_phi; | ||
|
||
// tout = t; | ||
// ptout = pt; | ||
|
||
// assign updated values | ||
x = xout; | ||
y = yout; | ||
t = tout; | ||
px = pxout; | ||
py = pyout; | ||
pt = ptout; | ||
|
||
// undo shift due to alignment errors of the element | ||
shift_out(x, y, px, py); | ||
} | ||
|
||
/** This pushes the reference particle. */ | ||
using Thin::operator(); | ||
|
||
amrex::ParticleReal m_phi; //! angle of rotation (rad) | ||
}; | ||
|
||
} // namespace impactx | ||
|
||
#endif // IMPACTX_PLANE_XYROT_H |
Oops, something went wrong.
Oops, something went wrong.
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
Check notice
Code scanning / CodeQL
Commented-out code Note