You can not select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
246 lines
8.9 KiB
246 lines
8.9 KiB
2 years ago
|
"""Function to get the default dictionary for fixtures."""
|
||
|
|
||
|
import argparse
|
||
|
import pathlib
|
||
|
|
||
|
import numpy
|
||
|
import shapely.geometry
|
||
|
import json
|
||
|
|
||
|
DEFAULT_TIMESTEP = 1e-2
|
||
|
DEFAULT_INITIAL_EPSILON = 1e-1
|
||
|
DEFAULT_RESTITUTION_COEFFICIENT = -1
|
||
|
DEFAULT_GRAVITY = [0.0, 0.0, 0.0]
|
||
|
DEFAULT_NUM_STEPS = 1000
|
||
|
DEFAULT_BARRIER_SOLVER_EB = 1e-6
|
||
|
DEFAULT_BARRIER_SOLVER_C = 0.1
|
||
|
DEFAULT_BARRIER_SOLVER_TINIT = 100
|
||
|
DEFAULT_BARRIER_SOLVER_TINC = 100
|
||
|
|
||
|
|
||
|
def generate_default_fixture() -> dict:
|
||
|
"""Create the default fixture as a dictionary."""
|
||
|
return {
|
||
|
"scene_type": "distance_barrier_rb_problem",
|
||
|
"max_iterations": DEFAULT_NUM_STEPS,
|
||
|
"timestep": DEFAULT_TIMESTEP,
|
||
|
"distance_barrier_constraint": {
|
||
|
"initial_barrier_activation_distance": DEFAULT_INITIAL_EPSILON,
|
||
|
"detection_method": "hash_grid"
|
||
|
},
|
||
|
"barrier_solver": {
|
||
|
"e_b": DEFAULT_BARRIER_SOLVER_EB,
|
||
|
"m": 1,
|
||
|
"t_init": DEFAULT_BARRIER_SOLVER_TINIT,
|
||
|
"t_inc": DEFAULT_BARRIER_SOLVER_TINC,
|
||
|
"c": DEFAULT_BARRIER_SOLVER_C,
|
||
|
"inner_solver": "newton_solver"
|
||
|
},
|
||
|
"rigid_body_problem": {
|
||
|
"gravity": DEFAULT_GRAVITY,
|
||
|
"coefficient_restitution": DEFAULT_RESTITUTION_COEFFICIENT,
|
||
|
"rigid_bodies": []
|
||
|
}
|
||
|
}
|
||
|
|
||
|
|
||
|
def generate_custom_fixture(args: argparse.Namespace) -> dict:
|
||
|
fixture = generate_default_fixture()
|
||
|
fixture["timestep"] = args.timestep
|
||
|
fixture["max_iterations"] = args.num_steps
|
||
|
fixture["distance_barrier_constraint"]["initial_barrier_activation_distance"] = (
|
||
|
args.init_epsilon)
|
||
|
fixture["rigid_body_problem"]["gravity"] = (args.gravity + 3 * [0])[:3]
|
||
|
fixture["rigid_body_problem"]["coefficient_restitution"] = (
|
||
|
args.restitution_coeff)
|
||
|
fixture["barrier_solver"]["e_b"] = args.eb
|
||
|
fixture["barrier_solver"]["c"] = args.c
|
||
|
fixture["barrier_solver"]["t_init"] = args.tinit
|
||
|
fixture["barrier_solver"]["t_inc"] = args.tinc
|
||
|
return fixture
|
||
|
|
||
|
|
||
|
def create_argument_parser(
|
||
|
description: str,
|
||
|
default_timestep: float = DEFAULT_TIMESTEP,
|
||
|
default_initial_epsilon: float = DEFAULT_INITIAL_EPSILON,
|
||
|
default_restitution_coefficient:
|
||
|
float = DEFAULT_RESTITUTION_COEFFICIENT,
|
||
|
default_gravity: list = DEFAULT_GRAVITY,
|
||
|
default_num_steps: int = DEFAULT_NUM_STEPS,
|
||
|
default_barrier_solver_eb: float = DEFAULT_BARRIER_SOLVER_EB,
|
||
|
default_barrier_solver_c: float = DEFAULT_BARRIER_SOLVER_C,
|
||
|
default_barrier_solver_tinit: float = DEFAULT_BARRIER_SOLVER_TINIT,
|
||
|
default_barrier_solver_tinc: float = DEFAULT_BARRIER_SOLVER_TINC,
|
||
|
) -> argparse.ArgumentParser:
|
||
|
parser = argparse.ArgumentParser(description=description)
|
||
|
parser.add_argument("--time-step",
|
||
|
type=float,
|
||
|
default=default_timestep,
|
||
|
dest="timestep",
|
||
|
help="length of the time-step (Δt)")
|
||
|
parser.add_argument("--num-steps",
|
||
|
type=int,
|
||
|
default=default_num_steps,
|
||
|
help="number of time-steps to take")
|
||
|
parser.add_argument("--init-epsilon",
|
||
|
type=float,
|
||
|
default=default_initial_epsilon,
|
||
|
help="initial d̂ for the barrier")
|
||
|
parser.add_argument("--gravity",
|
||
|
type=float,
|
||
|
nargs=2,
|
||
|
default=default_gravity,
|
||
|
help="[x, y] vector for gravitational acceleration")
|
||
|
parser.add_argument("--cor",
|
||
|
"--restitution-coefficient",
|
||
|
type=float,
|
||
|
dest="restitution_coeff",
|
||
|
default=default_restitution_coefficient,
|
||
|
help="coefficient of restitution")
|
||
|
parser.add_argument("--out-path",
|
||
|
metavar="path/to/output.json",
|
||
|
type=pathlib.Path,
|
||
|
default=None,
|
||
|
help="path to save the fixture")
|
||
|
parser.add_argument("--eb", type=float,
|
||
|
default=default_barrier_solver_eb)
|
||
|
parser.add_argument("--c", type=float,
|
||
|
default=default_barrier_solver_c)
|
||
|
parser.add_argument("--tinit", type=float,
|
||
|
default=default_barrier_solver_tinit)
|
||
|
parser.add_argument("--tinc", type=float,
|
||
|
default=default_barrier_solver_tinc)
|
||
|
return parser
|
||
|
|
||
|
|
||
|
def generate_regular_ngon_vertices(n: int, radius: float) -> numpy.ndarray:
|
||
|
"""Generate the vertices of a regular N-gon centered at zero."""
|
||
|
angles = (numpy.arange(n, dtype=float) * 2 * numpy.pi / n).reshape(-1, 1)
|
||
|
return radius * numpy.hstack([numpy.cos(angles), numpy.sin(angles)])
|
||
|
|
||
|
|
||
|
def generate_ngon_edges(n: int) -> numpy.ndarray:
|
||
|
"""Generate the edges of a N-gon."""
|
||
|
indices = numpy.arange(n).reshape(-1, 1)
|
||
|
return numpy.hstack([indices, numpy.roll(indices, -1)])
|
||
|
|
||
|
|
||
|
def generate_rectangle_vertices(hx: float, hy: float, center: numpy.ndarray,
|
||
|
angle: float) -> numpy.ndarray:
|
||
|
"""Generate a rectangle polygon."""
|
||
|
points = numpy.array([[hx, hy], [-hx, hy], [-hx, -hy], [hx, -hy]])
|
||
|
points = points @ create_2D_rotation_matrix(angle).T
|
||
|
points += center
|
||
|
return points
|
||
|
|
||
|
|
||
|
def generate_rectangle(hx: float, hy: float, center: numpy.ndarray,
|
||
|
angle: float) -> shapely.geometry.Polygon:
|
||
|
"""Generate a rectangle polygon."""
|
||
|
return shapely.geometry.Polygon(
|
||
|
generate_rectangle_vertices(hx, hy, center, angle))
|
||
|
|
||
|
|
||
|
def create_2D_rotation_matrix(theta: float) -> numpy.ndarray:
|
||
|
c, s = numpy.cos(theta), numpy.sin(theta)
|
||
|
return numpy.array([[c, -s], [s, c]])
|
||
|
|
||
|
|
||
|
def generate_walls_body(hx: float, hy: float, center: numpy.ndarray,
|
||
|
thickness: float) -> dict:
|
||
|
"""Generate a rigid body dictionary for walls."""
|
||
|
assert (thickness > 0)
|
||
|
# Inner vertices of the wall
|
||
|
inner_vertices = generate_rectangle_vertices(hx, hy, center, 0)
|
||
|
|
||
|
# Outer vertices of the wall
|
||
|
diag_thickness = thickness / numpy.sin(numpy.pi / 4)
|
||
|
outer_vertices = inner_vertices + thickness * numpy.array(
|
||
|
[[1, 1], [-1, 1], [-1, -1], [1, -1]])
|
||
|
|
||
|
# Combined vertices (inner vertices should be CW)
|
||
|
vertices = numpy.append(inner_vertices[::-1], outer_vertices, axis=0)
|
||
|
|
||
|
# Polygons of the wall
|
||
|
polygons = numpy.array([[
|
||
|
outer_vertices[i], outer_vertices[(i + 1) % 4],
|
||
|
inner_vertices[(i + 1) % 4], inner_vertices[i]
|
||
|
] for i in range(4)])
|
||
|
# Check that the polygons are all oriented counter-clockwise
|
||
|
for polygon in polygons:
|
||
|
assert is_polygon_ccw(polygon)
|
||
|
|
||
|
quad_edges = generate_ngon_edges(4)
|
||
|
edges = numpy.append(quad_edges, quad_edges + 4, axis=0)
|
||
|
return {
|
||
|
"vertices": vertices.tolist(),
|
||
|
"polygons": polygons.tolist(),
|
||
|
"edges": edges.tolist(),
|
||
|
"oriented": True,
|
||
|
"velocity": [0.0, 0.0, 0.0],
|
||
|
"is_dof_fixed": [True, True, True]
|
||
|
}
|
||
|
|
||
|
|
||
|
def print_args(args: argparse.Namespace) -> None:
|
||
|
""" Print the arguments."""
|
||
|
for k, v in args.__dict__.items():
|
||
|
print(f"{k}: {v}")
|
||
|
|
||
|
|
||
|
def save_fixture(fixture, fname):
|
||
|
with open(fname, "w") as outfile:
|
||
|
json.dump(fixture, outfile, indent=None, separators=(',', ':'))
|
||
|
|
||
|
|
||
|
def get_fixture_dir_path() -> pathlib.Path:
|
||
|
return pathlib.Path(__file__).resolve().parents[2] / "fixtures"
|
||
|
|
||
|
|
||
|
def get_meshes_dir_path() -> pathlib.Path:
|
||
|
return pathlib.Path(__file__).resolve().parents[2] / "meshes"
|
||
|
|
||
|
|
||
|
def is_polygon_ccw(vertices):
|
||
|
"""
|
||
|
Check if a polygon's vertices are in counter-clockwise order.
|
||
|
|
||
|
Uses the Shoelace Formula:
|
||
|
https://bit.ly/2t78ytv
|
||
|
https://en.wikipedia.org/wiki/Shoelace_formula
|
||
|
"""
|
||
|
winding_number = 0
|
||
|
for i in range(vertices.shape[0]):
|
||
|
x2_m_x1 = vertices[(i + 1) % vertices.shape[0], 0] - vertices[i, 0]
|
||
|
y2_p_y1 = vertices[(i + 1) % vertices.shape[0], 1] + vertices[i, 1]
|
||
|
winding_number += x2_m_x1 * y2_p_y1
|
||
|
return winding_number < 0
|
||
|
|
||
|
|
||
|
def compute_regular_ngon_area(vertices):
|
||
|
n = vertices.shape[0]
|
||
|
assert n > 3
|
||
|
side_length = numpy.linalg.norm(vertices[1] - vertices[0])
|
||
|
return n * side_length**2 / (4 * numpy.tan(numpy.pi / n)) # cm²
|
||
|
|
||
|
|
||
|
def generate_box_body(hx: float, hy: float, center: list, angle: float,
|
||
|
mass: float) -> dict:
|
||
|
vertices = generate_rectangle_vertices(hx, hy, [0, 0], 0)
|
||
|
edges = generate_ngon_edges(4)
|
||
|
area = 4 * hx * hy # m^2
|
||
|
density = mass / area # m^2
|
||
|
return {
|
||
|
"vertices": vertices.tolist(),
|
||
|
"polygons": [vertices.tolist()],
|
||
|
"edges": edges.tolist(),
|
||
|
"density": density,
|
||
|
"is_dof_fixed": numpy.zeros(3, dtype=bool).tolist(),
|
||
|
"oriented": True,
|
||
|
"linear_velocity": numpy.zeros(2).tolist(),
|
||
|
"angular_velocity": numpy.zeros(1).tolist(),
|
||
|
"position": center,
|
||
|
"rotation": [angle]
|
||
|
}
|