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.

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

feat: Elastic energy #10

Open
wants to merge 12 commits into
base: conditional_kernel
Choose a base branch
from
Open

Conversation

Vollkornaffe
Copy link
Collaborator

@Vollkornaffe Vollkornaffe commented Dec 6, 2023

This adds functions to compute the elastic energy density function from a given deformation gradient.

The new code is different from the existing one.
See pos_energy:

fn pos_energy(&self, particle: &Particle) -> Real {
self.pos_energy(particle.deformation_gradient, particle.elastic_hardening)
}

We may want to refactor.

@Vollkornaffe Vollkornaffe self-assigned this Dec 6, 2023
@Vollkornaffe Vollkornaffe changed the base branch from master to conditional_kernel December 6, 2023 18:06
@@ -17,6 +17,21 @@ impl ConstitutiveModel for CorotatedLinearElasticity {
)
}

// https://www.math.ucla.edu/~cffjiang/research/mpmcourse/mpmcourse.pdf#subsection.6.3
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
// https://www.math.ucla.edu/~cffjiang/research/mpmcourse/mpmcourse.pdf#subsection.6.3
// https://www.math.ucla.edu/~cffjiang/research/mpmcourse/mpmcourse.pdf#subsection.6.3
// This is essentially the same calculation as `[pos_energy]` but without ignoring the negative part.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I added a reference to the pos_energy instead.

Comment on lines 21 to 32
fn elastic_energy_density(&self, deformation_gradient: Matrix<Real>) -> Real {
let singular_values = deformation_gradient
.svd_unordered(false, false)
.singular_values;
let determinant: Real = singular_values.iter().product();

self.mu
* singular_values
.iter()
.map(|sigma| (sigma - 1.).powi(2))
.sum::<Real>()
+ self.lambda / 2. * (determinant - 1.).powi(2)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I would recommend taking the whole &Particle as argument, for consistency with pos_energy and to take hardening into account the same way as in pos_energy.

Can be in another PR.

Comment on lines 21 to 28
fn elastic_energy_density(&self, deformation_gradient: Matrix<Real>) -> Real {
let determinant_log = deformation_gradient.determinant().ln();
self.mu / 2.
* ((deformation_gradient.transpose() * deformation_gradient).trace() - DIM as Real)
- self.mu * determinant_log
+ self.lambda / 2. * determinant_log.powi(2)
}

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Same remark regarding taking hardening into account.

@@ -76,7 +76,7 @@ impl CorotatedLinearElasticity {

pub fn pos_energy(&self, deformation_gradient: Matrix<Real>, elastic_hardening: Real) -> Real {
let j = deformation_gradient.determinant();
let mut pos_def = deformation_gradient.svd_unordered(true, true);
let mut pos_def = deformation_gradient.svd_unordered(true, true); // TODO: why compute U and V?
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Probably just a refactoring leftover. It can be false, false.

// General elastic density function: https://dl.acm.org/doi/pdf/10.1145/3306346.3322949#section.5 (8) and (9)
// With degradation: https://dl.acm.org/doi/pdf/10.1145/3306346.3322949#subsection.3.2 (3)
// With hardening: https://www.math.ucla.edu/~cffjiang/research/mpmcourse/mpmcourse.pdf#subsection.6.5 (87)
#[allow(non_snake_case)]
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Let’s not allow non_snake_case. It messes up with syntax highlighting, and is very unconventional in rust, even though I agree it makes it closer to the mathematical notation.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We're now using more descriptive names (with proper snake case) while keeping aliases matching the paper in the comments.

Comment on lines 105 to 107
let Psi_mu =
|F: Matrix<Real>| hardened_mu / 2. * ((F.transpose() * F).trace() - DIM as Real);
let Psi_mu = Psi_mu(J_alpha * F);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The closure makes the code more difficult to read, especially since it’s only called once and bound to a variable with the same name.

@Vollkornaffe
Copy link
Collaborator Author

Once we've converged on style for the neo-hookean one, I'll update the corotated code to match.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

2 participants