Add trailing commas; clean up formatting
All checks were successful
/ test (pull_request) Successful in 3m55s

This commit is contained in:
Aaron Fenyes 2025-08-01 13:50:51 -07:00
parent 2eba80fb69
commit af59166906
12 changed files with 235 additions and 213 deletions

View file

@ -16,7 +16,7 @@ pub fn sphere(center_x: f64, center_y: f64, center_z: f64, radius: f64) -> DVect
center_y / radius,
center_z / radius,
0.5 / radius,
0.5 * (center_norm_sq / radius - radius)
0.5 * (center_norm_sq / radius - radius),
])
}
@ -30,7 +30,7 @@ pub fn sphere_with_offset(dir_x: f64, dir_y: f64, dir_z: f64, off: f64, curv: f6
norm_sp * dir_y,
norm_sp * dir_z,
0.5 * curv,
off * (1.0 + 0.5 * off * curv)
off * (1.0 + 0.5 * off * curv),
])
}
@ -53,7 +53,7 @@ pub fn project_point_to_normalized(rep: &mut DVector<f64>) {
pub struct MatrixEntry {
index: (usize, usize),
value: f64
value: f64,
}
pub struct PartialMatrix(Vec<MatrixEntry>);
@ -65,7 +65,7 @@ impl PartialMatrix {
pub fn push(&mut self, row: usize, col: usize, value: f64) {
let PartialMatrix(entries) = self;
entries.push(MatrixEntry { index: (row, col), value: value });
entries.push(MatrixEntry { index: (row, col), value });
}
pub fn push_sym(&mut self, row: usize, col: usize, value: f64) {
@ -135,22 +135,26 @@ impl<'a> IntoIterator for &'a PartialMatrix {
pub struct ConfigSubspace {
assembly_dim: usize,
basis_std: Vec<DMatrix<f64>>,
basis_proj: Vec<DMatrix<f64>>
basis_proj: Vec<DMatrix<f64>>,
}
impl ConfigSubspace {
pub fn zero(assembly_dim: usize) -> ConfigSubspace {
ConfigSubspace {
assembly_dim: assembly_dim,
assembly_dim,
basis_proj: Vec::new(),
basis_std: Vec::new()
basis_std: Vec::new(),
}
}
// approximate the kernel of a symmetric endomorphism of the configuration
// space for `assembly_dim` elements. we consider an eigenvector to be part
// of the kernel if its eigenvalue is smaller than the constant `THRESHOLD`
fn symmetric_kernel(a: DMatrix<f64>, proj_to_std: DMatrix<f64>, assembly_dim: usize) -> ConfigSubspace {
fn symmetric_kernel(
a: DMatrix<f64>,
proj_to_std: DMatrix<f64>,
assembly_dim: usize,
) -> ConfigSubspace {
// find a basis for the kernel. the basis is expressed in the projection
// coordinates, and it's orthonormal with respect to the projection
// inner product
@ -170,7 +174,7 @@ impl ConfigSubspace {
const ELEMENT_DIM: usize = 5;
const UNIFORM_DIM: usize = 4;
ConfigSubspace {
assembly_dim: assembly_dim,
assembly_dim,
basis_std: basis_std.column_iter().map(
|v| Into::<DMatrix<f64>>::into(
v.reshape_generic(Dyn(ELEMENT_DIM), Dyn(assembly_dim))
@ -180,7 +184,7 @@ impl ConfigSubspace {
|v| Into::<DMatrix<f64>>::into(
v.reshape_generic(Dyn(UNIFORM_DIM), Dyn(assembly_dim))
)
).collect()
).collect(),
}
}
@ -214,9 +218,9 @@ pub struct DescentHistory {
pub config: Vec<DMatrix<f64>>,
pub scaled_loss: Vec<f64>,
pub neg_grad: Vec<DMatrix<f64>>,
pub hess_eigvals: Vec::<DVector<f64>>,
pub hess_eigvals: Vec<DVector<f64>>,
pub base_step: Vec<DMatrix<f64>>,
pub backoff_steps: Vec<i32>
pub backoff_steps: Vec<i32>,
}
impl DescentHistory {
@ -246,7 +250,7 @@ impl ConstraintProblem {
ConstraintProblem {
gram: PartialMatrix::new(),
frozen: PartialMatrix::new(),
guess: DMatrix::<f64>::zeros(ELEMENT_DIM, element_count)
guess: DMatrix::<f64>::zeros(ELEMENT_DIM, element_count),
}
}
@ -255,7 +259,7 @@ impl ConstraintProblem {
ConstraintProblem {
gram: PartialMatrix::new(),
frozen: PartialMatrix::new(),
guess: DMatrix::from_columns(guess_columns)
guess: DMatrix::from_columns(guess_columns),
}
}
}
@ -269,25 +273,21 @@ lazy_static! {
0.0, 1.0, 0.0, 0.0, 0.0,
0.0, 0.0, 1.0, 0.0, 0.0,
0.0, 0.0, 0.0, 0.0, -2.0,
0.0, 0.0, 0.0, -2.0, 0.0
0.0, 0.0, 0.0, -2.0, 0.0,
]);
}
struct SearchState {
config: DMatrix<f64>,
err_proj: DMatrix<f64>,
loss: f64
loss: f64,
}
impl SearchState {
fn from_config(gram: &PartialMatrix, config: DMatrix<f64>) -> SearchState {
let err_proj = gram.sub_proj(&(config.tr_mul(&*Q) * &config));
let loss = err_proj.norm_squared();
SearchState {
config: config,
err_proj: err_proj,
loss: loss
}
SearchState { config, err_proj, loss }
}
}
@ -314,7 +314,7 @@ pub fn local_unif_to_std(v: DVectorView<f64>) -> DMatrix<f64> {
curv, 0.0, 0.0, 0.0, v[0],
0.0, curv, 0.0, 0.0, v[1],
0.0, 0.0, curv, 0.0, v[2],
0.0, 0.0, 0.0, 0.0, 1.0
0.0, 0.0, 0.0, 0.0, 1.0,
])
} else {
// `v` represents a sphere. the normalization condition says that the
@ -323,7 +323,7 @@ pub fn local_unif_to_std(v: DVectorView<f64>) -> DMatrix<f64> {
curv, 0.0, 0.0, 0.0, v[0],
0.0, curv, 0.0, 0.0, v[1],
0.0, 0.0, curv, 0.0, v[2],
curv*v[0], curv*v[1], curv*v[2], curv*v[3], curv*v[4] + 1.0
curv*v[0], curv*v[1], curv*v[2], curv*v[3], curv*v[4] + 1.0,
])
}
}
@ -336,7 +336,7 @@ fn seek_better_config(
base_target_improvement: f64,
min_efficiency: f64,
backoff: f64,
max_backoff_steps: i32
max_backoff_steps: i32,
) -> Option<(SearchState, i32)> {
let mut rate = 1.0;
for backoff_steps in 0..max_backoff_steps {
@ -354,12 +354,12 @@ fn seek_better_config(
// a first-order neighborhood of a configuration
pub struct ConfigNeighborhood {
pub config: DMatrix<f64>,
pub nbhd: ConfigSubspace
pub nbhd: ConfigSubspace,
}
pub struct Realization {
pub result: Result<ConfigNeighborhood, String>,
pub history: DescentHistory
pub history: DescentHistory,
}
// seek a matrix `config` that matches the partial matrix `problem.frozen` and
@ -373,12 +373,10 @@ pub fn realize_gram(
backoff: f64,
reg_scale: f64,
max_descent_steps: i32,
max_backoff_steps: i32
max_backoff_steps: i32,
) -> Realization {
// destructure the problem data
let ConstraintProblem {
gram, guess, frozen
} = problem;
let ConstraintProblem { gram, guess, frozen } = problem;
// start the descent history
let mut history = DescentHistory::new();
@ -391,10 +389,10 @@ pub fn realize_gram(
let result = Ok(
ConfigNeighborhood {
config: guess.clone(),
nbhd: ConfigSubspace::zero(0)
nbhd: ConfigSubspace::zero(0),
}
);
return Realization { result, history }
return Realization { result, history };
}
// find the dimension of the search space
@ -475,7 +473,7 @@ pub fn realize_gram(
Some(cholesky) => cholesky,
None => return Realization {
result: Err("Cholesky decomposition failed".to_string()),
history
history,
}
};
let base_step_stacked = hess_cholesky.solve(&neg_grad_stacked);
@ -485,16 +483,16 @@ pub fn realize_gram(
// use backtracking line search to find a better configuration
if let Some((better_state, backoff_steps)) = seek_better_config(
gram, &state, &base_step, neg_grad.dot(&base_step),
min_efficiency, backoff, max_backoff_steps
min_efficiency, backoff, max_backoff_steps,
) {
state = better_state;
history.backoff_steps.push(backoff_steps);
} else {
return Realization {
result: Err("Line search failed".to_string()),
history
}
};
history,
};
}
}
let result = if state.loss < tol {
// express the uniform basis in the standard basis
@ -539,7 +537,7 @@ pub mod examples {
[
sphere(0.0, 0.0, 0.0, 15.0),
sphere(0.0, 0.0, -9.0, 5.0),
sphere(0.0, 0.0, 11.0, 3.0)
sphere(0.0, 0.0, 11.0, 3.0),
].into_iter().chain(
(1..=6).map(
|k| {
@ -598,7 +596,7 @@ pub mod examples {
point(0.0, 0.0, 0.0),
point(ang_hor.cos(), ang_hor.sin(), 0.0),
point(x_vert, y_vert, -0.5),
point(x_vert, y_vert, 0.5)
point(x_vert, y_vert, 0.5),
]
}
).collect::<Vec<_>>().as_slice()
@ -641,15 +639,15 @@ mod tests {
MatrixEntry { index: (0, 0), value: 14.0 },
MatrixEntry { index: (0, 2), value: 28.0 },
MatrixEntry { index: (1, 1), value: 42.0 },
MatrixEntry { index: (1, 2), value: 49.0 }
MatrixEntry { index: (1, 2), value: 49.0 },
]);
let config = DMatrix::<f64>::from_row_slice(2, 3, &[
1.0, 2.0, 3.0,
4.0, 5.0, 6.0
4.0, 5.0, 6.0,
]);
let expected_result = DMatrix::<f64>::from_row_slice(2, 3, &[
14.0, 2.0, 28.0,
4.0, 42.0, 49.0
4.0, 42.0, 49.0,
]);
assert_eq!(frozen.freeze(&config), expected_result);
}
@ -660,15 +658,15 @@ mod tests {
MatrixEntry { index: (0, 0), value: 19.0 },
MatrixEntry { index: (0, 2), value: 39.0 },
MatrixEntry { index: (1, 1), value: 59.0 },
MatrixEntry { index: (1, 2), value: 69.0 }
MatrixEntry { index: (1, 2), value: 69.0 },
]);
let attempt = DMatrix::<f64>::from_row_slice(2, 3, &[
1.0, 2.0, 3.0,
4.0, 5.0, 6.0
4.0, 5.0, 6.0,
]);
let expected_result = DMatrix::<f64>::from_row_slice(2, 3, &[
18.0, 0.0, 36.0,
0.0, 54.0, 63.0
0.0, 54.0, 63.0,
]);
assert_eq!(target.sub_proj(&attempt), expected_result);
}
@ -686,7 +684,7 @@ mod tests {
DMatrix::from_columns(&[
sphere(1.0, 0.0, 0.0, a),
sphere(-0.5, a, 0.0, a),
sphere(-0.5, -a, 0.0, a)
sphere(-0.5, -a, 0.0, a),
])
};
let state = SearchState::from_config(&gram, config);
@ -700,7 +698,7 @@ mod tests {
fn frozen_entry_test() {
let mut problem = ConstraintProblem::from_guess(&[
point(0.0, 0.0, 2.0),
sphere(0.0, 0.0, 0.0, 0.95)
sphere(0.0, 0.0, 0.0, 0.95),
]);
for j in 0..2 {
for k in j..2 {
@ -744,7 +742,7 @@ mod tests {
let mut problem = ConstraintProblem::from_guess(&[
sphere(0.0, 0.0, 0.0, -2.0),
sphere(0.0, 0.0, 1.0, 1.0),
sphere(0.0, 0.0, -1.0, 1.0)
sphere(0.0, 0.0, -1.0, 1.0),
]);
for j in 0..3 {
for k in j..3 {
@ -774,8 +772,8 @@ mod tests {
DMatrix::<f64>::from_column_slice(UNIFORM_DIM, assembly_dim, &[
0.0, 0.0, 0.0, 0.0,
0.0, 0.0, -0.5, -0.5,
0.0, 0.0, -0.5, 0.5
])
0.0, 0.0, -0.5, 0.5,
]),
];
let tangent_motions_std = vec![
basis_matrix((0, 1), element_dim, assembly_dim),
@ -785,8 +783,8 @@ mod tests {
DMatrix::<f64>::from_column_slice(element_dim, assembly_dim, &[
0.0, 0.0, 0.0, 0.00, 0.0,
0.0, 0.0, -1.0, -0.25, -1.0,
0.0, 0.0, -1.0, 0.25, 1.0
])
0.0, 0.0, -1.0, 0.25, 1.0,
]),
];
// confirm that the dimension of the tangent space is no greater than
@ -862,10 +860,10 @@ mod tests {
DVector::from_column_slice(&[0.0, 0.0, 5.0, 0.0]),
DVector::from_column_slice(&[0.0, 0.0, 1.0, 0.0]),
DVector::from_column_slice(&[-vel_vert_x, -vel_vert_y, -3.0, 0.0]),
DVector::from_column_slice(&[vel_vert_x, vel_vert_y, -3.0, 0.0])
DVector::from_column_slice(&[vel_vert_x, vel_vert_y, -3.0, 0.0]),
]
}
).collect::<Vec<_>>()
).collect::<Vec<_>>(),
];
let tangent_motions_std = tangent_motions_unif.iter().map(
|motion| DMatrix::from_columns(
@ -898,7 +896,7 @@ mod tests {
0.0, 1.0, 0.0, 0.0, dis[1],
0.0, 0.0, 1.0, 0.0, dis[2],
2.0*dis[0], 2.0*dis[1], 2.0*dis[2], 1.0, dis.norm_squared(),
0.0, 0.0, 0.0, 0.0, 1.0
0.0, 0.0, 0.0, 0.0, 1.0,
])
}
@ -910,7 +908,7 @@ mod tests {
const SCALED_TOL: f64 = 1.0e-12;
let mut problem_orig = ConstraintProblem::from_guess(&[
sphere(0.0, 0.0, 0.5, 1.0),
sphere(0.0, 0.0, -0.5, 1.0)
sphere(0.0, 0.0, -0.5, 1.0),
]);
problem_orig.gram.push_sym(0, 0, 1.0);
problem_orig.gram.push_sym(1, 1, 1.0);
@ -928,13 +926,13 @@ mod tests {
let a = 0.5 * FRAC_1_SQRT_2;
DMatrix::from_columns(&[
sphere(a, 0.0, 7.0 + a, 1.0),
sphere(-a, 0.0, 7.0 - a, 1.0)
sphere(-a, 0.0, 7.0 - a, 1.0),
])
};
let problem_tfm = ConstraintProblem {
gram: problem_orig.gram,
frozen: problem_orig.frozen,
guess: guess_tfm,
frozen: problem_orig.frozen
};
let Realization { result: result_tfm, history: history_tfm } = realize_gram(
&problem_tfm, SCALED_TOL, 0.5, 0.9, 1.1, 200, 110
@ -962,7 +960,7 @@ mod tests {
0.0, 1.0, 0.0, 0.0, 0.0,
FRAC_1_SQRT_2, 0.0, FRAC_1_SQRT_2, 0.0, 0.0,
0.0, 0.0, 0.0, 1.0, 0.0,
0.0, 0.0, 0.0, 0.0, 1.0
0.0, 0.0, 0.0, 0.0, 1.0,
]);
let transl = translation(Vector3::new(0.0, 0.0, 7.0));
let motion_proj_tfm = transl * rot * motion_orig_proj;