Potentials¶
Barrier Potential¶
- class ipctk.BarrierPotential¶
Bases:
pybind11_object
Public Methods:
__init__
(*args, **kwargs)Overloaded function.
__call__
(*args, **kwargs)Overloaded function.
gradient
(*args, **kwargs)Overloaded function.
hessian
(*args, **kwargs)Overloaded function.
shape_derivative
(*args, **kwargs)Overloaded function.
Inherited from
pybind11_object
-
__annotations__ =
{}
¶
- __call__(*args, **kwargs)¶
Overloaded function.
__call__(self: ipctk.BarrierPotential, collisions: ipctk.Collisions, mesh: ipc::CollisionMesh, vertices: numpy.ndarray[numpy.float64[m, n]]) -> float
Compute the barrier potential for a set of collisions.
- Parameters:
collisions: The set of collisions. mesh: The collision mesh. vertices: Vertices of the collision mesh.
- Returns:
The sum of all barrier potentials (not scaled by the barrier stiffness).
__call__(self: ipctk.BarrierPotential, collision: ipctk.Collision, x: numpy.ndarray[numpy.float64[m, 1]]) -> float
Compute the potential for a single collision.
- Parameters:
collision: The collision. x: The collision stencil’s degrees of freedom.
- Returns:
The potential.
- __init__(*args, **kwargs)¶
Overloaded function.
__init__(self: ipctk.BarrierPotential, dhat: float, use_physical_barrier: bool = False) -> None
Construct a barrier potential.
- Parameters:
dhat: The activation distance of the barrier.
__init__(self: ipctk.BarrierPotential, barrier: ipctk.Barrier, dhat: float, use_physical_barrier: bool = False) -> None
Construct a barrier potential.
- Parameters:
barrier: The barrier function. dhat: The activation distance of the barrier.
-
__module__ =
'ipctk'
¶
- property barrier : ipctk.Barrier¶
Barrier activation distance.
- property dhat : float¶
Barrier activation distance.
- gradient(*args, **kwargs)¶
Overloaded function.
gradient(self: ipctk.BarrierPotential, collisions: ipctk.Collisions, mesh: ipc::CollisionMesh, vertices: numpy.ndarray[numpy.float64[m, n]]) -> numpy.ndarray[numpy.float64[m, 1]]
Compute the gradient of the barrier potential.
- Parameters:
collisions: The set of collisions. mesh: The collision mesh. vertices: Vertices of the collision mesh.
- Returns:
The gradient of all barrier potentials (not scaled by the barrier stiffness). This will have a size of |vertices|.
gradient(self: ipctk.BarrierPotential, collision: ipctk.Collision, x: numpy.ndarray[numpy.float64[m, 1]]) -> numpy.ndarray[numpy.float64[m, 1]]
Compute the gradient of the potential for a single collision.
- Parameters:
collision: The collision. x: The collision stencil’s degrees of freedom.
- Returns:
The gradient of the potential.
- hessian(*args, **kwargs)¶
Overloaded function.
hessian(self: ipctk.BarrierPotential, collisions: ipctk.Collisions, mesh: ipc::CollisionMesh, vertices: numpy.ndarray[numpy.float64[m, n]], project_hessian_to_psd: ipctk.PSDProjectionMethod = <PSDProjectionMethod.NONE: 0>) -> scipy.sparse.csc_matrix[numpy.float64]
Compute the hessian of the barrier potential.
- Parameters:
collisions: The set of collisions. mesh: The collision mesh. vertices: Vertices of the collision mesh. project_hessian_to_psd: Make sure the hessian is positive semi-definite.
- Returns:
The hessian of all barrier potentials (not scaled by the barrier stiffness). This will have a size of |vertices|x|vertices|.
hessian(self: ipctk.BarrierPotential, collision: ipctk.Collision, x: numpy.ndarray[numpy.float64[m, 1]], project_hessian_to_psd: ipctk.PSDProjectionMethod = <PSDProjectionMethod.NONE: 0>) -> numpy.ndarray[numpy.float64[m, n]]
Compute the hessian of the potential for a single collision.
- Parameters:
collision: The collision. x: The collision stencil’s degrees of freedom.
- Returns:
The hessian of the potential.
- shape_derivative(*args, **kwargs)¶
Overloaded function.
shape_derivative(self: ipctk.BarrierPotential, collisions: ipctk.Collisions, mesh: ipc::CollisionMesh, vertices: numpy.ndarray[numpy.float64[m, n]]) -> scipy.sparse.csc_matrix[numpy.float64]
Compute the shape derivative of the potential.
std::runtime_error If the collision collisions were not built with shape derivatives enabled.
- Parameters:
collisions: The set of collisions. mesh: The collision mesh. vertices: Vertices of the collision mesh.
- Returns:
The derivative of the force with respect to X, the rest vertices.
shape_derivative(self: ipc::DistanceBasedPotential, collision: ipctk.Collision, vertex_ids: Annotated[list[int], FixedSize(4)], rest_positions: numpy.ndarray[numpy.float64[m, 1]], positions: numpy.ndarray[numpy.float64[m, 1]]) -> scipy.sparse.csc_matrix[numpy.float64]
Compute the shape derivative of the potential for a single collision.
- Parameters:
collision: The collision. vertex_ids: The collision stencil’s vertex ids. rest_positions: The collision stencil’s rest positions. positions: The collision stencil’s positions. ,out]: out Store the triplets of the shape derivative here.
-
__annotations__ =
Friction Potential¶
- class ipctk.FrictionPotential¶
Bases:
pybind11_object
Public Data Attributes:
The smooth friction mollifier parameter \(\epsilon_{v}\).
Public Methods:
__init__
(self, epsv)Construct a friction potential.
__call__
(*args, **kwargs)Overloaded function.
gradient
(*args, **kwargs)Overloaded function.
hessian
(*args, **kwargs)Overloaded function.
force
(*args, **kwargs)Overloaded function.
force_jacobian
(*args, **kwargs)Overloaded function.
Inherited from
pybind11_object
- class DiffWRT¶
Bases:
pybind11_object
Members:
REST_POSITIONS : Differentiate w.r.t. rest positions
LAGGED_DISPLACEMENTS : Differentiate w.r.t. lagged displacements
VELOCITIES : Differentiate w.r.t. current velocities
-
LAGGED_DISPLACEMENTS =
<DiffWRT.LAGGED_DISPLACEMENTS: 1>
¶
-
REST_POSITIONS =
<DiffWRT.REST_POSITIONS: 0>
¶
-
VELOCITIES =
<DiffWRT.VELOCITIES: 2>
¶
-
__annotations__ =
{}
¶
- __eq__(self, other: object) bool ¶
- __getstate__(self) int ¶
- __hash__(self) int ¶
- __index__(self) int ¶
- __init__(self, value: int)¶
- __int__(self) int ¶
-
__members__ =
{'LAGGED_DISPLACEMENTS': <DiffWRT.LAGGED_DISPLACEMENTS: 1>, 'REST_POSITIONS': <DiffWRT.REST_POSITIONS: 0>, 'VELOCITIES': <DiffWRT.VELOCITIES: 2>}
¶
-
__module__ =
'ipctk'
¶
- __ne__(self, other: object) bool ¶
- __repr__(self) str ¶
- __setstate__(self, state: int) None ¶
- __str__(self) str ¶
- property name : str¶
- property value : int¶
-
LAGGED_DISPLACEMENTS =
-
LAGGED_DISPLACEMENTS =
<DiffWRT.LAGGED_DISPLACEMENTS: 1>
¶
-
REST_POSITIONS =
<DiffWRT.REST_POSITIONS: 0>
¶
-
VELOCITIES =
<DiffWRT.VELOCITIES: 2>
¶
-
__annotations__ =
{}
¶
- __call__(*args, **kwargs)¶
Overloaded function.
__call__(self: ipctk.FrictionPotential, collisions: ipctk.FrictionCollisions, mesh: ipc::CollisionMesh, velocities: numpy.ndarray[numpy.float64[m, n]]) -> float
Compute the friction dissipative potential for a set of collisions.
- Parameters:
collisions: The set of collisions. mesh: The collision mesh. velocities: Velocities of the collision mesh.
- Returns:
The sum of all friction dissipative potentials.
__call__(self: ipctk.FrictionPotential, collision: ipctk.FrictionCollision, x: numpy.ndarray[numpy.float64[m, 1]]) -> float
Compute the potential for a single collision.
- Parameters:
collision: The collision. x: The collision stencil’s degrees of freedom.
- Returns:
The potential.
-
__module__ =
'ipctk'
¶
- property epsv : float¶
The smooth friction mollifier parameter \(\epsilon_{v}\).
- force(*args, **kwargs)¶
Overloaded function.
force(self: ipctk.FrictionPotential, collisions: ipctk.FrictionCollisions, mesh: ipc::CollisionMesh, rest_positions: numpy.ndarray[numpy.float64[m, n]], lagged_displacements: numpy.ndarray[numpy.float64[m, n]], velocities: numpy.ndarray[numpy.float64[m, n]], barrier_potential: ipctk.BarrierPotential, barrier_stiffness: float, dmin: float = 0, no_mu: bool = False) -> numpy.ndarray[numpy.float64[m, 1]]
Compute the friction force for all collisions.
- Parameters:
collisions: The set of collisions. mesh: The collision mesh. rest_positions: Rest positions of the vertices (rowwise). lagged_displacements: Previous displacements of the vertices (rowwise). velocities: Current velocities of the vertices (rowwise). barrier_potential: Barrier potential (used for normal force magnitude). barrier_stiffness: Barrier stiffness (used for normal force magnitude). dmin: Minimum distance (used for normal force magnitude). no_mu: whether to not multiply by mu
- Returns:
The friction force.
force(self: ipctk.FrictionPotential, collision: ipctk.FrictionCollision, rest_positions: numpy.ndarray[numpy.float64[m, 1]], lagged_displacements: numpy.ndarray[numpy.float64[m, 1]], velocities: numpy.ndarray[numpy.float64[m, 1]], barrier_potential: ipctk.BarrierPotential, barrier_stiffness: float, dmin: float = 0, no_mu: bool = False) -> numpy.ndarray[numpy.float64[m, 1]]
Compute the friction force for a single collision.
- Parameters:
collision: The collision rest_positions: Rest positions of the vertices (rowwise). lagged_displacements: Previous displacements of the vertices (rowwise). velocities: Current displacements of the vertices (rowwise). barrier_potential: Barrier potential (used for normal force magnitude). barrier_stiffness: Barrier stiffness (used for normal force magnitude). dmin: Minimum distance (used for normal force magnitude). no_mu: Whether to not multiply by mu
- Returns:
Friction force
- force_jacobian(*args, **kwargs)¶
Overloaded function.
force_jacobian(self: ipctk.FrictionPotential, collisions: ipctk.FrictionCollisions, mesh: ipc::CollisionMesh, rest_positions: numpy.ndarray[numpy.float64[m, n]], lagged_displacements: numpy.ndarray[numpy.float64[m, n]], velocities: numpy.ndarray[numpy.float64[m, n]], barrier_potential: ipctk.BarrierPotential, barrier_stiffness: float, wrt: ipctk.FrictionPotential.DiffWRT, dmin: float = 0) -> scipy.sparse.csc_matrix[numpy.float64]
Compute the Jacobian of the friction force for all collisions.
- Parameters:
collisions: The set of collisions. mesh: The collision mesh. rest_positions: Rest positions of the vertices (rowwise). lagged_displacements: Previous displacements of the vertices (rowwise). velocities: Current displacements of the vertices (rowwise). barrier_potential: Barrier potential (used for normal force magnitude). barrier_stiffness: Barrier stiffness (used for normal force magnitude). wrt: The variable to take the derivative with respect to. dmin: Minimum distance (used for normal force magnitude).
- Returns:
The Jacobian of the friction force wrt the velocities.
force_jacobian(self: ipctk.FrictionPotential, collision: ipctk.FrictionCollision, rest_positions: numpy.ndarray[numpy.float64[m, 1]], lagged_displacements: numpy.ndarray[numpy.float64[m, 1]], velocities: numpy.ndarray[numpy.float64[m, 1]], barrier_potential: ipctk.BarrierPotential, barrier_stiffness: float, wrt: ipctk.FrictionPotential.DiffWRT, dmin: float = 0) -> numpy.ndarray[numpy.float64[m, n]]
Compute the friction force Jacobian.
- Parameters:
collision: The collision rest_positions: Rest positions of the vertices (rowwise). lagged_displacements: Previous displacements of the vertices (rowwise). velocities: Current displacements of the vertices (rowwise). barrier_potential: Barrier potential (used for normal force magnitude). barrier_stiffness: Barrier stiffness (used for normal force magnitude). wrt: Variable to differentiate the friction force with respect to. dmin: Minimum distance (used for normal force magnitude).
- Returns:
Friction force Jacobian
- gradient(*args, **kwargs)¶
Overloaded function.
gradient(self: ipctk.FrictionPotential, collisions: ipctk.FrictionCollisions, mesh: ipc::CollisionMesh, velocities: numpy.ndarray[numpy.float64[m, n]]) -> numpy.ndarray[numpy.float64[m, 1]]
Compute the gradient of the friction dissipative potential.
- Parameters:
collisions: The set of collisions. mesh: The collision mesh. velocities: Velocities of the collision mesh.
- Returns:
The gradient of all friction dissipative potentials. This will have a size of |velocities|.
gradient(self: ipctk.FrictionPotential, collision: ipctk.FrictionCollision, x: numpy.ndarray[numpy.float64[m, 1]]) -> numpy.ndarray[numpy.float64[m, 1]]
Compute the gradient of the potential for a single collision.
- Parameters:
collision: The collision. x: The collision stencil’s degrees of freedom.
- Returns:
The gradient of the potential.
- hessian(*args, **kwargs)¶
Overloaded function.
hessian(self: ipctk.FrictionPotential, collisions: ipctk.FrictionCollisions, mesh: ipc::CollisionMesh, velocities: numpy.ndarray[numpy.float64[m, n]], project_hessian_to_psd: ipctk.PSDProjectionMethod = <PSDProjectionMethod.NONE: 0>) -> scipy.sparse.csc_matrix[numpy.float64]
Compute the hessian of the friction dissipative potential.
- Parameters:
collisions: The set of collisions. mesh: The collision mesh. velocities: Velocities of the collision mesh. project_hessian_to_psd: Make sure the hessian is positive semi-definite.
- Returns:
The hessian of all friction dissipative potentials. This will have a size of |velocities|×|velocities|.
hessian(self: ipctk.FrictionPotential, collision: ipctk.FrictionCollision, x: numpy.ndarray[numpy.float64[m, 1]], project_hessian_to_psd: ipctk.PSDProjectionMethod = <PSDProjectionMethod.NONE: 0>) -> numpy.ndarray[numpy.float64[m, n]]
Compute the hessian of the potential for a single collision.
- Parameters:
collision: The collision. x: The collision stencil’s degrees of freedom.
- Returns:
The hessian of the potential.