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.

Sign up for GitHub

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

Move force details into models #172

Open
wants to merge 2 commits into
base: main
Choose a base branch
from
Open

Move force details into models #172

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
e5973cd
Remove unused default constructors and set_param
streeve Jul 11, 2024
22ce383
Add model interfaces for force/energy calculations
streeve Jan 9, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
15 changes: 0 additions & 15 deletions src/CabanaPD_ForceModels.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -21,12 +21,9 @@ struct BaseForceModel
{
double delta;

BaseForceModel(){};
BaseForceModel( const double _delta )
: delta( _delta ){};

BaseForceModel( const BaseForceModel& model ) { delta = model.delta; }

// No-op for temperature.
KOKKOS_INLINE_FUNCTION
void thermalStretch( double&, const int, const int ) const {}
Expand Down Expand Up @@ -55,12 +52,6 @@ struct BaseTemperatureModel
temperature = model.temperature;
}

void set_param( const double _alpha, const double _temp0 )
{
alpha = _alpha;
temp0 = _temp0;
}

void update( const TemperatureType _temp ) { temperature = _temp; }

// Update stretch with temperature effects.
Expand All @@ -87,12 +78,6 @@ struct BaseDynamicTemperatureModel
BaseDynamicTemperatureModel( const double _delta, const double _kappa,
const double _cp,
const bool _constant_microconductivity = true )
{
set_param( _delta, _kappa, _cp, _constant_microconductivity );
}

void set_param( const double _delta, const double _kappa, const double _cp,
const bool _constant_microconductivity )
{
delta = _delta;
kappa = _kappa;
Expand Down
1 change: 0 additions & 1 deletion src/CabanaPD_Solver.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -464,7 +464,6 @@ class SolverNoFracture
// Optional modules.
std::shared_ptr<heat_transfer_type> heat_transfer;
std::shared_ptr<contact_type> contact;
contact_model_type contact_model;

// Output files.
std::string output_file;
Expand Down
18 changes: 9 additions & 9 deletions src/force/CabanaPD_ContactModels.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,6 @@ struct ContactModel
double delta;
double Rc;

ContactModel(){};
// PD horizon
// Contact radius
ContactModel( const double _delta, const double _Rc )
Expand All @@ -51,22 +50,23 @@ struct NormalRepulsionModel : public ContactModel
double c;
double K;

NormalRepulsionModel(){};
NormalRepulsionModel( const double delta, const double Rc, const double _K )
: ContactModel( delta, Rc )
, K( _K )
{
set_param( delta, Rc, K );
}

void set_param( const double _delta, const double _Rc, const double _K )
{
delta = _delta;
Rc = _Rc;
K = _K;
// This could inherit from PMB (same c)
c = 18.0 * K / ( pi * delta * delta * delta * delta );
}

KOKKOS_INLINE_FUNCTION
auto forceCoeff( const double r, const double vol ) const
{
// Contact "stretch"
const double sc = ( r - Rc ) / delta;
// Normal repulsion uses a 15 factor compared to the PMB force
return 15.0 * c * sc * vol;
}
};

template <>
Expand Down
58 changes: 40 additions & 18 deletions src/force/CabanaPD_ForceModels_LPS.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -40,27 +40,56 @@ struct ForceModel<LPS, Elastic, NoFracture> : public BaseForceModel
const int _influence = 0 )
: base_type( _delta )
, influence_type( _influence )
, K( _K )
, G( _G )
{
set_param( _delta, _K, _G );
}

void set_param( const double _delta, const double _K, const double _G )
{
delta = _delta;
K = _K;
G = _G;

theta_coeff = 3.0 * K - 5.0 * G;
s_coeff = 15.0 * G;
}

KOKKOS_INLINE_FUNCTION double influence_function( double xi ) const
KOKKOS_INLINE_FUNCTION double influenceFunction( double xi ) const
{
if ( influence_type == 1 )
return 1.0 / xi;
else
return 1.0;
}

KOKKOS_INLINE_FUNCTION auto weightedVolume( const double xi,
const double vol ) const
{
return influenceFunction( xi ) * xi * xi * vol;
}

KOKKOS_INLINE_FUNCTION auto dilatation( const double s, const double xi,
const double vol,
const double m_i ) const
{
double theta_i = influenceFunction( xi ) * s * xi * xi * vol;
return 3.0 * theta_i / m_i;
}

KOKKOS_INLINE_FUNCTION auto forceCoeff( const double s, const double xi,
const double vol, const double m_i,
const double m_j,
const double theta_i,
const double theta_j ) const
{
return ( theta_coeff * ( theta_i / m_i + theta_j / m_j ) +
s_coeff * s * ( 1.0 / m_i + 1.0 / m_j ) ) *
influenceFunction( xi ) * xi * vol;
}

KOKKOS_INLINE_FUNCTION
auto energy( const double s, const double xi, const double vol,
const double m_i, const double theta_i,
const double num_bonds ) const
{
return 1.0 / num_bonds * 0.5 * theta_coeff / 3.0 *
( theta_i * theta_i ) +
0.5 * ( s_coeff / m_i ) * influenceFunction( xi ) * s * s * xi *
xi * vol;
}
};

template <>
Expand All @@ -85,15 +114,8 @@ struct ForceModel<LPS, Elastic, Fracture>
ForceModel( const double _delta, const double _K, const double _G,
const double _G0, const int _influence = 0 )
: base_type( _delta, _K, _G, _influence )
, G0( _G0 )
{
set_param( _delta, _K, _G, _G0 );
}

void set_param( const double _delta, const double _K, const double _G,
const double _G0 )
{
base_type::set_param( _delta, _K, _G );
G0 = _G0;
if ( influence_type == 1 )
{
s0 = Kokkos::sqrt( 5.0 * G0 / 9.0 / K / delta ); // 1/xi
Expand Down
81 changes: 15 additions & 66 deletions src/force/CabanaPD_ForceModels_PMB.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -34,24 +34,25 @@ struct ForceModel<PMB, Elastic, NoFracture, TemperatureIndependent>
double c;
double K;

ForceModel( const double delta, const double K )
ForceModel( const double delta, const double _K )
: base_type( delta )
, K( _K )
{
set_param( delta, K );
c = 18.0 * K / ( pi * delta * delta * delta * delta );
}

ForceModel( const ForceModel& model )
: base_type( model )
KOKKOS_INLINE_FUNCTION
auto forceCoeff( const double s, const double vol ) const
{
c = model.c;
K = model.K;
return c * s * vol;
}

void set_param( const double _delta, const double _K )
KOKKOS_INLINE_FUNCTION
auto energy( const double s, const double xi, const double vol ) const
{
delta = _delta;
K = _K;
c = 18.0 * K / ( pi * delta * delta * delta * delta );
// 0.25 factor is due to 1/2 from outside the integral and 1/2 from
// the integrand (pairwise potential).
return 0.25 * c * s * s * xi * vol;
}
};

Expand All @@ -71,24 +72,10 @@ struct ForceModel<PMB, Elastic, Fracture, TemperatureIndependent>
double s0;
double bond_break_coeff;

ForceModel( const double delta, const double K, const double G0 )
ForceModel( const double delta, const double K, const double _G0 )
: base_type( delta, K )
, G0( _G0 )
{
set_param( delta, K, G0 );
}

ForceModel( const ForceModel& model )
: base_type( model )
{
G0 = model.G0;
s0 = model.s0;
bond_break_coeff = model.bond_break_coeff;
}

void set_param( const double _delta, const double _K, const double _G0 )
{
base_type::set_param( _delta, _K );
G0 = _G0;
s0 = Kokkos::sqrt( 5.0 * G0 / 9.0 / K / delta );
bond_break_coeff = ( 1.0 + s0 ) * ( 1.0 + s0 );
}
Expand Down Expand Up @@ -167,14 +154,6 @@ struct ForceModel<PMB, Elastic, NoFracture, TemperatureDependent,
: base_type( _delta, _K )
, base_temperature_type( _temp, _alpha, _temp0 )
{
set_param( _delta, _K, _alpha, _temp0 );
}

void set_param( const double _delta, const double _K, const double _alpha,
const double _temp0 )
{
base_type::set_param( _delta, _K );
base_temperature_type::set_param( _alpha, _temp0 );
}
};

Expand Down Expand Up @@ -232,14 +211,6 @@ struct ForceModel<PMB, Elastic, Fracture, TemperatureDependent, TemperatureType>
: base_type( _delta, _K, _G0 )
, base_temperature_type( _temp, _alpha, _temp0 )
{
set_param( _delta, _K, _G0, _alpha, _temp0 );
}

void set_param( const double _delta, const double _K, const double _G0,
const double _alpha, const double _temp0 )
{
base_type::set_param( _delta, _K, _G0 );
base_temperature_type::set_param( _alpha, _temp0 );
}

KOKKOS_INLINE_FUNCTION
Expand Down Expand Up @@ -301,17 +272,6 @@ struct ForceModel<PMB, Elastic, NoFracture, DynamicTemperature, TemperatureType>
, base_temperature_type( _delta, _kappa, _cp,
_constant_microconductivity )
{
set_param( _delta, _K, _kappa, _cp, _alpha, _temp0,
_constant_microconductivity );
}

void set_param( const double _delta, const double _K, const double _kappa,
const double _cp, const double _alpha, const double _temp0,
const bool _constant_microconductivity )
{
base_type::set_param( _delta, _K, _alpha, _temp0 );
base_temperature_type::set_param( _delta, _kappa, _cp,
_constant_microconductivity );
}
};

Expand Down Expand Up @@ -364,20 +324,9 @@ struct ForceModel<PMB, Fracture, DynamicTemperature, TemperatureType>
const double _temp0 = 0.0,
const bool _constant_microconductivity = true )
: base_type( _delta, _K, _G0, _temp, _alpha, _temp0 )
, base_temperature_type( _delta, _kappa, _cp,
_constant_microconductivity )
{
set_param( _delta, _K, _G0, _kappa, _cp, _alpha, _temp0 );
base_temperature_type::set_param( _delta, _kappa, _cp,
_constant_microconductivity );
}

void set_param( const double _delta, const double _K, const double _G0,
const double _kappa, const double _cp, const double _alpha,
const double _temp0,
const bool _constant_microconductivity )
{
base_type::set_param( _delta, _K, _G0, _alpha, _temp0 );
base_temperature_type::set_param( _delta, _kappa, _cp,
_constant_microconductivity );
}
};

Expand Down
13 changes: 4 additions & 9 deletions src/force/CabanaPD_Force_Contact.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -69,16 +69,15 @@ class Force<MemorySpace, NormalRepulsionModel>
const ParticleType& particles,
ParallelType& neigh_op_tag )
{
auto delta = _model.delta;
auto Rc = _model.Rc;
auto c = _model.c;
auto model = _model;
const auto vol = particles.sliceVolume();
const auto y = particles.sliceCurrentPosition();
const int n_frozen = particles.frozenOffset();
const int n_local = particles.localOffset();

_neigh_timer.start();
_neigh_list.build( y, n_frozen, n_local, Rc, 1.0, mesh_min, mesh_max );
_neigh_list.build( y, n_frozen, n_local, model.Rc, 1.0, mesh_min,
mesh_max );
_neigh_timer.stop();

auto contact_full = KOKKOS_LAMBDA( const int i, const int j )
Expand All @@ -91,11 +90,7 @@ class Force<MemorySpace, NormalRepulsionModel>
double rx, ry, rz;
getDistanceComponents( x, u, i, j, xi, r, s, rx, ry, rz );

// Contact "stretch"
const double sc = ( r - Rc ) / delta;

// Normal repulsion uses a 15 factor compared to the PMB force
const double coeff = 15 * c * sc * vol( j );
const double coeff = model.forceCoeff( r, vol( j ) );
fcx_i = coeff * rx / r;
fcy_i = coeff * ry / r;
fcz_i = coeff * rz / r;
Expand Down
Loading
Loading