Skip to content
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.

Sign up for GitHub

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

Jump to bottom

機能の分離 #6

Merged
merged 1 commit into from
Aug 8, 2023
Merged

機能の分離 #6

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
86 changes: 56 additions & 30 deletions src/kd_tree.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,17 +1,52 @@
use super::gen_grid2d;
#[derive(Debug, Clone, PartialEq, PartialOrd)]
pub struct Grid2D {
pub x: f64,
pub y: f64,
}

#[derive(Debug, Clone)]
pub struct Points2D {
pub points: Vec<Grid2D>,
}

#[allow(dead_code)]
#[derive(Debug, Clone)]
pub struct KDTree {
id: usize,
position: gen_grid2d::Grid2D,
position: Grid2D,
left: Option<Box<KDTree>>,
right: Option<Box<KDTree>>,
}

impl Grid2D {
#[allow(dead_code)]
pub fn new(x_: f64, y_: f64) -> Self {
Grid2D { x: x_, y: y_ }
}

#[allow(dead_code)]
pub fn distance_square(&self, point: &Grid2D) -> f64 {
let dx2 = (self.x - point.x) * (self.x - point.x);
let dy2 = (self.y - point.y) * (self.y - point.y);
dx2 + dy2
}
}

impl Points2D {
#[allow(dead_code)]
pub fn new() -> Self {
Points2D { points: vec![] }
}

#[allow(dead_code)]
pub fn push(&mut self, x_r: f64, y_r: f64) {
self.points.push(Grid2D { x: x_r, y: y_r });
}
}

impl KDTree {
#[allow(dead_code)]
fn new(vector: &gen_grid2d::Grid2D, id_: usize) -> Self {
fn new(vector: &Grid2D, id_: usize) -> Self {
Self {
id: id_,
position: vector.clone(),
Expand Down Expand Up @@ -43,20 +78,14 @@ impl KDTree {
}

#[allow(dead_code)]
pub fn neighbor_search(&self, x: &gen_grid2d::Grid2D, radius: f64) -> Vec<usize> {
pub fn neighbor_search(&self, x: &Grid2D, radius: f64) -> Vec<usize> {
let mut near = vec![0; 0];
self.search_points_id(x, radius, &mut near, 0);
near.clone()
}

#[allow(dead_code)]
fn search_points_id(
&self,
x: &gen_grid2d::Grid2D,
radius: f64,
near: &mut Vec<usize>,
mut depth: i32,
) {
fn search_points_id(&self, x: &Grid2D, radius: f64, near: &mut Vec<usize>, mut depth: i32) {
let axis = depth % 2;
let r_self = self.position.distance_square(x).sqrt();
if r_self < radius {
Expand Down Expand Up @@ -140,7 +169,7 @@ impl KDTree {
}

#[allow(dead_code)]
fn insert(&mut self, point: &gen_grid2d::Grid2D, mut depth: i32, id: usize) {
fn insert(&mut self, point: &Grid2D, mut depth: i32, id: usize) {
let axis = depth % 2;

match axis {
Expand Down Expand Up @@ -198,7 +227,7 @@ impl KDTree {
}

#[allow(dead_code)]
fn create_kd_tree(&mut self, vec: &gen_grid2d::Points2D) -> Self {
fn create_kd_tree(&mut self, vec: &Points2D) -> Self {
let depth = 0;
for i in 1..vec.points.len() {
self.insert(&vec.points[i], depth, i);
Expand All @@ -217,11 +246,8 @@ impl KDTree {
}

#[allow(dead_code)]
pub fn construct_kd_tree(vec: &gen_grid2d::Points2D) -> KDTree {
let mut tree = KDTree::new(
&gen_grid2d::Grid2D::new(vec.points[0].x, vec.points[0].y),
0,
);
pub fn construct_kd_tree(vec: &Points2D) -> KDTree {
let mut tree = KDTree::new(&Grid2D::new(vec.points[0].x, vec.points[0].y), 0);
tree.create_kd_tree(vec)
}

Expand All @@ -247,7 +273,7 @@ mod tests {
let mut rng: rand::rngs::StdRng = rand::SeedableRng::from_seed(seed);
let num_point: usize = 600;

let mut vec = gen_grid2d::Points2D::new();
let mut vec = Points2D::new();

for _ in 0..num_point {
let x_r = 2.0 * (rng.gen::<f64>() - 0.5);
Expand All @@ -268,7 +294,7 @@ mod tests {
let mut rng: rand::rngs::StdRng = rand::SeedableRng::from_seed(seed);
let num_point: usize = 10;

let mut vec = gen_grid2d::Points2D::new();
let mut vec = Points2D::new();

for _ in 0..num_point {
let x_r = 2.0 * (rng.gen::<f64>() - 0.5);
Expand All @@ -277,7 +303,7 @@ mod tests {
}

let tree = KDTree::construct_kd_tree(&mut vec);
let center = gen_grid2d::Grid2D { x: 0.0, y: 0.0 };
let center = Grid2D { x: 0.0, y: 0.0 };
let radius = 0.4;
let near = tree.neighbor_search(&center, radius);

Expand All @@ -291,7 +317,7 @@ mod tests {
let mut rng: rand::rngs::StdRng = rand::SeedableRng::from_seed(seed);
let num_point: usize = 10;

let mut vec = gen_grid2d::Points2D::new();
let mut vec = Points2D::new();

for _ in 0..num_point {
let x_r = 2.0 * (rng.gen::<f64>() - 0.5);
Expand All @@ -300,7 +326,7 @@ mod tests {
}

let tree = KDTree::construct_kd_tree(&mut vec);
let center = gen_grid2d::Grid2D { x: 0.4, y: 0.3 };
let center = Grid2D { x: 0.4, y: 0.3 };
let radius = 0.5;
let near = tree.neighbor_search(&center, radius);

Expand All @@ -317,7 +343,7 @@ mod tests {
let mut rng: rand::rngs::StdRng = rand::SeedableRng::from_seed(seed);
let num_point: usize = 10;

let mut vec = gen_grid2d::Points2D::new();
let mut vec = Points2D::new();

for _ in 0..num_point {
let x_r = 2.0 * (rng.gen::<f64>() - 0.5);
Expand All @@ -327,7 +353,7 @@ mod tests {

let tree = KDTree::construct_kd_tree(&mut vec);

let center = gen_grid2d::Grid2D { x: 0.4, y: 0.3 };
let center = Grid2D { x: 0.4, y: 0.3 };
let radius = 0.5;
let mut near = vec![0 as usize; 0];
tree.search_points_id(&center, radius, &mut near, 0);
Expand All @@ -342,7 +368,7 @@ mod tests {
let mut rnd = rand::thread_rng();
for _ in 0..10 {
let num_point: i32 = rnd.gen_range(0..100000);
let mut vec = gen_grid2d::Points2D::new();
let mut vec = Points2D::new();
for _ in 0..num_point {
let x_r = rnd.gen::<i32>() as f64;
let y_r = rnd.gen::<i32>() as f64;
Expand All @@ -355,13 +381,13 @@ mod tests {
}

#[allow(dead_code)]
fn benchmark_pre(size: i32) -> gen_grid2d::Points2D {
fn benchmark_pre(size: i32) -> Points2D {
use rand::Rng;

let mut rnd = rand::thread_rng();
let max = size;
let num_point: i32 = rnd.gen_range(0..max);
let mut vec = gen_grid2d::Points2D::new();
let mut vec = Points2D::new();
for _ in 0..num_point {
let x_r = rnd.gen::<i32>() as f64 / max as f64;
let y_r = rnd.gen::<i32>() as f64 / max as f64;
Expand All @@ -371,9 +397,9 @@ mod tests {
}

#[allow(dead_code)]
fn benchmark(vec: &gen_grid2d::Points2D) {
fn benchmark(vec: &Points2D) {
let tree = KDTree::construct_kd_tree(vec);
let center = gen_grid2d::Grid2D::new(0.0, 0.0);
let center = Grid2D::new(0.0, 0.0);
tree.neighbor_search(&center, 0.01);
}

Expand Down
120 changes: 38 additions & 82 deletions src/gen_grid2d.rs → src/lennard_jones_potential.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,45 +1,10 @@
use super::kd_tree;

#[derive(Debug, Clone, PartialEq, PartialOrd)]
pub struct Grid2D {
pub x: f64,
pub y: f64,
}

impl Grid2D {
#[allow(dead_code)]
pub fn new(x_: f64, y_: f64) -> Self {
Grid2D { x: x_, y: y_ }
}

#[allow(dead_code)]
pub fn distance_square(&self, point: &Grid2D) -> f64 {
let dx2 = (self.x - point.x) * (self.x - point.x);
let dy2 = (self.y - point.y) * (self.y - point.y);
dx2 + dy2
}
}

#[derive(Debug, Clone)]
pub struct Points2D {
pub points: Vec<Grid2D>,
}

impl Points2D {
#[allow(dead_code)]
pub fn new() -> Self {
Points2D { points: vec![] }
}

#[allow(dead_code)]
pub fn push(&mut self, x_r: f64, y_r: f64) {
self.points.push(Grid2D { x: x_r, y: y_r });
}

impl kd_tree::Points2D {
#[allow(dead_code)]
pub fn euler_step_by_near_points(
&mut self,
boundary: &Points2D,
boundary: &kd_tree::Points2D,
tree: &kd_tree::KDTree,
radius: f64,
) {
Expand Down Expand Up @@ -72,7 +37,7 @@ impl Points2D {
}

#[allow(dead_code)]
pub fn euler_step(&mut self, boundary: &Points2D) {
pub fn euler_step(&mut self, boundary: &kd_tree::Points2D) {
let dt = 1.0e-3;
for i in 0..self.points.len() {
let mut x = self.points[i].x - dt * self.lennard_jones_potential_deriv(i, boundary).x;
Expand All @@ -99,10 +64,10 @@ impl Points2D {
pub fn lennard_jones_potential_deriv_by_near_points(
&mut self,
index: usize,
boundary: &Points2D,
boundary: &kd_tree::Points2D,
tree: &kd_tree::KDTree,
radius: f64,
) -> Grid2D {
) -> kd_tree::Grid2D {
let epsilon: f64 = 1.0;
let sigma: f64 = 0.9 * (2.0_f64).powf(-1.0 / 6.0) / (self.points.len() as f64).sqrt();
let mut f_x = 0.0;
Expand All @@ -123,16 +88,8 @@ impl Points2D {
let r8 = r4 * r4;
let r13 = r8 * r4 * r;
let r7 = r2 * r4 * r;
f_x += 4.0
* epsilon
* (12.0 * sigma12 / r13
- 6.0 * sigma6 / r7)
* dx;
f_y += 4.0
* epsilon
* (12.0 * sigma12 / r13
- 6.0 * sigma6 / r7)
* dy;
f_x += 4.0 * epsilon * (12.0 * sigma12 / r13 - 6.0 * sigma6 / r7) * dx;
f_y += 4.0 * epsilon * (12.0 * sigma12 / r13 - 6.0 * sigma6 / r7) * dy;
}
}
for i in 0..boundary.points.len() {
Expand All @@ -144,54 +101,53 @@ impl Points2D {
let r8 = r4 * r4;
let r13 = r8 * r4 * r;
let r7 = r2 * r4 * r;
f_x += 4.0
* epsilon
* (12.0 * sigma12 / r13 - 6.0 * sigma6 / r7)
* dx;
f_y += 4.0
* epsilon
* (12.0 * sigma12 / r13 - 6.0 * sigma6 / r7)
* dy;
f_x += 4.0 * epsilon * (12.0 * sigma12 / r13 - 6.0 * sigma6 / r7) * dx;
f_y += 4.0 * epsilon * (12.0 * sigma12 / r13 - 6.0 * sigma6 / r7) * dy;
}
Grid2D { x: f_x, y: f_y }
kd_tree::Grid2D { x: f_x, y: f_y }
}

#[allow(dead_code)]
pub fn lennard_jones_potential_deriv(&mut self, index: usize, boundary: &Points2D) -> Grid2D {
pub fn lennard_jones_potential_deriv(
&mut self,
index: usize,
boundary: &kd_tree::Points2D,
) -> kd_tree::Grid2D {
let epsilon: f64 = 1.0;
let sigma: f64 = 0.9 * (2.0_f64).powf(-1.0 / 6.0) / (self.points.len() as f64).powf(0.5);
let sigma: f64 = 0.9 * (2.0_f64).powf(-1.0 / 6.0) / (self.points.len() as f64).sqrt();
let mut f_x = 0.0;
let mut f_y = 0.0;
let sigma2 = sigma * sigma;
let sigma4 = sigma2 * sigma2;
let sigma8 = sigma4 * sigma4;
let sigma6 = sigma2 * sigma4;
let sigma12 = sigma8 * sigma4;
for i in 0..self.points.len() {
if i != index {
let dx = self.points[index].x - self.points[i].x;
let dy = self.points[index].y - self.points[i].y;
let r = (dx.powf(2.0) + dy.powf(2.0)).powf(0.5);
f_x += 4.0
* epsilon
* (12.0 * sigma.powf(12.0) / r.powf(13.0)
- 6.0 * sigma.powf(6.0) / r.powf(7.0))
* dx;
f_y += 4.0
* epsilon
* (12.0 * sigma.powf(12.0) / r.powf(13.0)
- 6.0 * sigma.powf(6.0) / r.powf(7.0))
* dy;
let r = (dx * dx + dy * dy).sqrt();
let r2 = r * r;
let r4 = r2 * r2;
let r8 = r4 * r4;
let r13 = r8 * r4 * r;
let r7 = r2 * r4 * r;
f_x += 4.0 * epsilon * (12.0 * sigma12 / r13 - 6.0 * sigma6 / r7) * dx;
f_y += 4.0 * epsilon * (12.0 * sigma12 / r13 - 6.0 * sigma6 / r7) * dy;
}
}
for i in 0..boundary.points.len() {
let dx = self.points[index].x - boundary.points[i].x;
let dy = self.points[index].y - boundary.points[i].y;
let r = (dx.powf(2.0) + dy.powf(2.0)).powf(0.5);
f_x += 4.0
* epsilon
* (12.0 * sigma.powf(12.0) / r.powf(13.0) - 6.0 * sigma.powf(6.0) / r.powf(7.0))
* dx;
f_y += 4.0
* epsilon
* (12.0 * sigma.powf(12.0) / r.powf(13.0) - 6.0 * sigma.powf(6.0) / r.powf(7.0))
* dy;
let r = (dx * dx + dy * dy).sqrt();
let r2 = r * r;
let r4 = r2 * r2;
let r8 = r4 * r4;
let r13 = r8 * r4 * r;
let r7 = r2 * r4 * r;
f_x += 4.0 * epsilon * (12.0 * sigma12 / r13 - 6.0 * sigma6 / r7) * dx;
f_y += 4.0 * epsilon * (12.0 * sigma12 / r13 - 6.0 * sigma6 / r7) * dy;
}
Grid2D { x: f_x, y: f_y }
kd_tree::Grid2D { x: f_x, y: f_y }
}
}
Loading