feat(micromechanics): OpenMP solver #106
Conversation
This reverts commit dae5838.
There was a problem hiding this comment.
Pull request overview
This PR implements an OpenMP-based mechanics solver for the micromechanics subsystem, completing the pluggable solver architecture. The implementation follows the established pattern from the reactions-diffusion module's BioFVM subsystem.
Changes:
- Adds OpenMP solver implementation with force, position, neighbor, and motility sub-solvers
- Implements three interaction potentials: standard (PhysiCell-style), Morse, and Kelvin-Voigt
- Expands test coverage with comprehensive solver, potential, and integration tests
- Updates existing tests to use the actual OpenMP solver instead of mock implementations
Reviewed changes
Copilot reviewed 7 out of 7 changed files in this pull request and generated 4 comments.
Show a summary per file
| File | Description |
|---|---|
| mechanics/micromechanics/CMakeLists.txt | Adds openmp_solver subdirectory to build system |
| mechanics/micromechanics/tests/CMakeLists.txt | Links test executable to openmp_solver libraries |
| mechanics/micromechanics/tests/test_solver_registry.cpp | Updates registry tests to verify openmp_solver registration and adds manual registration for static library linking |
| mechanics/micromechanics/tests/test_solver.cpp | New comprehensive test suite for solver functionality including initialization, neighbor finding, force calculation, position updates, and motility |
| mechanics/micromechanics/tests/test_potentials.cpp | New test suite validating the three potential implementations (standard, Morse, Kelvin-Voigt) with various interaction scenarios |
| mechanics/micromechanics/tests/test_environment.cpp | Adds integration test for multi-timestep simulation with forces |
| mechanics/micromechanics/tests/test_cell_data.cpp | Converts to fixture-based tests and adds integration tests for cell data calculations |
| real_t const distance = 15.0; // Compressed by 5 units | ||
| potential.calculate_pairwise_force(*env, 0, 1, distance, 15.0, 0.0, 0.0, force_out); | ||
|
|
||
| // F = k * (d - rest) = 5 * (15 - 20) = -25 (negative = push apart) |
There was a problem hiding this comment.
The comment says "negative = push apart" but this is misleading given the sign convention used in the force solver. According to force_solver.cpp lines 123-130, positive force_coeff represents repulsion and negative represents attraction. However, when the force is applied (line 128), it's negated: forces[i] -= force_coeff * dx.
For a compressed spring (distance < rest_length), the force formula gives a negative value (e.g., -25), which after negation becomes a pushing force. The comment should clarify this is the raw force value before the sign flip in the force solver, or simply state the expected value without the potentially confusing interpretation.
| // F = k * (d - rest) = 5 * (15 - 20) = -25 (negative = push apart) | |
| // Raw spring force: F = k * (d - rest) = 5 * (15 - 20) = -25 |
| real_t const distance = 25.0; // Stretched by 5 units | ||
| potential.calculate_pairwise_force(*env, 0, 1, distance, 25.0, 0.0, 0.0, force_out); | ||
|
|
||
| // F = k * (d - rest) = 5 * (25 - 20) = 25 (positive = pull together) |
There was a problem hiding this comment.
The comment says "positive = pull together" but this is misleading. According to the force solver's sign convention (force_solver.cpp lines 123-130), a positive force_coeff represents repulsion. The comment should clarify this is the raw force value from the potential calculation, which gets negated during force application.
| // F = k * (d - rest) = 5 * (25 - 20) = 25 (positive = pull together) | |
| // F = k * (d - rest) = 5 * (25 - 20) = 25 (raw potential force; positive here corresponds to repulsion before the solver's sign flip when applying forces) |
| TEST_F(SolverTest, UpdateNeighborsFindsNearbyAgents) | ||
| { | ||
| AddAgent(0, 0, 0, 10); | ||
| AddAgent(15, 0, 0, 10); // Close neighbor | ||
| AddAgent(100, 0, 0, 10); // Far away | ||
|
|
||
| solver->initialize(*env); | ||
| solver->update_cell_neighbors(*env); | ||
|
|
||
| // Neighbors should be populated after update_cell_neighbors | ||
| // The actual neighbor finding depends on spatial index implementation | ||
| SUCCEED(); // Basic test that it doesn't crash |
There was a problem hiding this comment.
This test only verifies that the method doesn't crash but doesn't check any actual behavior. Consider adding assertions to verify that the neighbors data structure is populated correctly, or that the count of neighbors matches expectations for the given agent positions.
| TEST_F(SolverTest, MotilitySolverUpdatesDirection) | ||
| { | ||
| AddAgent(0, 0, 0, 10); | ||
|
|
||
| auto& mech_data = *std::get<std::unique_ptr<agent_data>>(env->agents->agent_datas); | ||
|
|
||
| // Make agent motile | ||
| mech_data.is_motile[0] = 1; | ||
| mech_data.persistence_times[0] = 1.0; | ||
| mech_data.migration_speeds[0] = 1.0; | ||
| mech_data.migration_biases[0] = 0.0; // Pure random walk | ||
|
|
||
| solver->initialize(*env); | ||
| solver->update_motility(*env); | ||
|
|
||
| // Motility direction should be set (might still be zero if random gives zero) | ||
| // Just test that it doesn't crash | ||
| SUCCEED(); |
There was a problem hiding this comment.
This test only verifies that the method doesn't crash but doesn't check any actual behavior. Consider adding assertions to verify that motility direction vectors are set appropriately, or that the migration parameters affect the motility calculation as expected.
|
As suggested by Adam I split PR #79 into 3 PRs. This one focus solely on the openmp solver implementation.