MDStressLab++
|
Computes deformation gradients and stress tensors using Moving Least Squares (MLS) interpolation on a spatial grid. More...
#include <Mls.h>
Public Member Functions | |
Mls (const BoxConfiguration &body, const Grid< Reference > *pgrid, double radiusMls, const std::string name) | |
Constructs the deformation gradient field from a body configuration and grid. | |
~Mls () | |
void | pushToCauchy (const std::vector< Matrix3d > &piolaStress, std::vector< Matrix3d > &cauchyStress) |
Converts Piola-Kirchhoff stress to Cauchy stress at grid points. | |
void | writeDeformationGradient () |
Writes the deformation gradient at each grid point to a file. | |
void | writeGridPushed () |
Writes the pushed (deformed) grid point coordinates to a file. | |
void | writePushedCauchyStress (std::vector< Matrix3d > &cauchyStress) |
Writes Cauchy stresses at grid points to a file. | |
Public Attributes | |
std::string | name |
Identifier name used for output file prefixes. | |
double | radiusMls |
Radius of influence for the MLS weighting function. | |
std::vector< Matrix3d > | deformationGradient |
Deformation gradient at each grid point (one per grid coordinate). | |
std::vector< Vector3d > | gridPushed |
Positions of grid points after applying MLS displacements. | |
Computes deformation gradients and stress tensors using Moving Least Squares (MLS) interpolation on a spatial grid.
This class implements a 3D MLS approach to compute the local deformation gradient and map Piola-Kirchhoff stress tensors to Cauchy stresses on a grid. It accounts for periodic boundary conditions and constructs local neighbor lists for accuracy.
Typical use: compute deformations and push stresses from atomistic data onto grid points.
Mls::Mls | ( | const BoxConfiguration & | body, |
const Grid< Reference > * | pgrid, | ||
double | radiusMls, | ||
const std::string | name | ||
) |
Constructs the deformation gradient field from a body configuration and grid.
Performs MLS interpolation to compute deformation gradients at each grid point using displacement fields derived from a reference and current configuration.
Handles periodic boundary conditions, builds neighbor lists, checks sanity of interpolation points, and stores the final deformation gradient and displaced grid points.
body | The configuration containing reference and current atom positions. |
pgrid | Pointer to the spatial grid (reference coordinates). |
radiusMls | The cutoff radius for MLS interpolation. |
name | An identifier used in output filenames. |
void Mls::pushToCauchy | ( | const std::vector< Matrix3d > & | piolaStress, |
std::vector< Matrix3d > & | cauchyStress | ||
) |
Converts Piola-Kirchhoff stress to Cauchy stress at grid points.
Uses the computed deformation gradients and their determinants to push stress tensors from the reference configuration to the current spatial grid.
piolaStress | Input Piola-Kirchhoff stresses at reference points. |
cauchyStress | Output Cauchy stresses on the deformed grid. |
void Mls::writeDeformationGradient | ( | ) |
Writes the deformation gradient at each grid point to a file.
File name is constructed as <name>.DeformationGradient
. Includes both grid positions and corresponding 3x3 deformation gradient matrices.
void Mls::writeGridPushed | ( | ) |
Writes the pushed (deformed) grid point coordinates to a file.
File name is <name>.pushedGrids
. One line per grid point with x/y/z coordinates.
void Mls::writePushedCauchyStress | ( | std::vector< Matrix3d > & | cauchyStress | ) |
Writes Cauchy stresses at grid points to a file.
File name is <name>.stress
. Outputs 6 unique components per symmetric 3x3 stress tensor.
cauchyStress | Vector of Cauchy stress tensors (must match grid size). |
std::vector<Matrix3d> Mls::deformationGradient |
std::vector<Vector3d> Mls::gridPushed |
std::string Mls::name |
double Mls::radiusMls |