Skip to content

Commit

Permalink
Revise
Browse files Browse the repository at this point in the history
  • Loading branch information
@TLCFEM
TLCFEM committed Sep 5, 2023
1 parent bfa1143 commit 2c1b789
Show file tree
Hide file tree
Showing 3 changed files with 56 additions and 48 deletions.
31 changes: 13 additions & 18 deletions Example/Material/ConcreteK4.supan
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,45 +14,40 @@ fix2 2 2 2 3
hdf5recorder 1 Element E 1
hdf5recorder 2 Element S 1

# step static 1
# set ini_step_size 1E-2
# set fixed_step_size 1
#
# displacement 1 0 -.04 1 3

step static 2
step static 1
set ini_step_size 1E-2
set fixed_step_size 1

displacement 2 0 .002 1 3
displacement 1 0 .002 1 3
converger RelIncreDisp 1 1E-8 10 1

step static 3
set ini_step_size 1E-3
step static 2
set ini_step_size 1E-1
set fixed_step_size 1

displacement 3 0 -.025 1 3
displacement 2 0 -.025 1 3

step static 4
set ini_step_size 1E-3
step static 3
set ini_step_size 1E-2
set fixed_step_size 1

displacement 4 0 .01 1 3
displacement 3 0 .01 1 3

step static 5
step static 4
set ini_step_size 1E-2
set fixed_step_size 1

displacement 5 0 -.02 1 3
displacement 4 0 -.02 1 3

analyze

# Node 3:
# Coordinate:
# 6.0000e+00 0.0000e+00
# Displacement:
# -1.0000e-01 0.0000e+00
# -3.3000e-02 0.0000e+00
# Resistance:
# -6.9916e+02 0.0000e+00
# -8.9176e+01 0.0000e+00
peek node 3

peek element 1 2
Expand Down
71 changes: 42 additions & 29 deletions Material/Material1D/Concrete/NonlinearK4.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -26,18 +26,22 @@ int NonlinearK4::compute_tension_branch() {

auto counter = 0u;
while(true) {
if(max_iteration == ++counter) return SUANPAN_FAIL;
if(max_iteration == ++counter) {
suanpan_error("Cannot converge within {} iterations.\n", max_iteration);
return SUANPAN_FAIL;
}

const auto backbone = compute_tension_backbone(kt);
const auto residual = fabs(trial_stress(0)) - backbone(0);

if(1u == counter && residual <= 0.) {
if(no_damage) return SUANPAN_SUCCESS;

const auto damage = compute_tension_damage(kt);
if(apply_damage) {
const auto damage = compute_tension_damage(kt);
const auto damage_factor = 1. - damage(0);

trial_stress *= 1. - damage(0);
trial_stiffness *= 1. - damage(0);
trial_stress *= damage_factor;
trial_stiffness *= damage_factor;
}

return SUANPAN_SUCCESS;
}
Expand All @@ -47,16 +51,18 @@ int NonlinearK4::compute_tension_branch() {
const auto error = fabs(incre);
suanpan_debug("Local tension iteration error: {:.5E}.\n", error);
if(error < tolerance || fabs(residual) < tolerance) {
trial_stiffness -= elastic_modulus / jacobian * elastic_modulus;

if(no_damage) return SUANPAN_SUCCESS;
const auto dgamma = elastic_modulus / jacobian;
trial_stiffness -= dgamma * elastic_modulus;

const auto damage = compute_tension_damage(kt);
if(apply_damage) {
const auto damage = compute_tension_damage(kt);
const auto damage_factor = 1. - damage(0);

trial_stiffness *= 1. - damage(0);
trial_stiffness -= trial_stress * damage(1) * elastic_modulus / jacobian;
trial_stiffness *= damage_factor;
trial_stiffness -= trial_stress * damage(1) * dgamma;

trial_stress *= 1. - damage(0);
trial_stress *= damage_factor;
}

return SUANPAN_SUCCESS;
}
Expand All @@ -77,18 +83,22 @@ int NonlinearK4::compute_compression_branch() {

auto counter = 0u;
while(true) {
if(max_iteration == ++counter) return SUANPAN_FAIL;
if(max_iteration == ++counter) {
suanpan_error("Cannot converge within {} iterations.\n", max_iteration);
return SUANPAN_FAIL;
}

const auto backbone = compute_compression_backbone(kc);
const auto residual = fabs(trial_stress(0)) - backbone(0);

if(1u == counter && residual <= 0.) {
if(no_damage) return SUANPAN_SUCCESS;

const auto damage = compute_compression_damage(kc);
if(apply_damage) {
const auto damage = compute_compression_damage(kc);
const auto damage_factor = 1. - damage(0);

trial_stress *= 1. - damage(0);
trial_stiffness *= 1. - damage(0);
trial_stress *= damage_factor;
trial_stiffness *= damage_factor;
}

return SUANPAN_SUCCESS;
}
Expand All @@ -98,16 +108,18 @@ int NonlinearK4::compute_compression_branch() {
const auto error = fabs(incre);
suanpan_debug("Local compression iteration error: {:.5E}.\n", error);
if(error < tolerance || fabs(residual) < tolerance) {
trial_stiffness -= elastic_modulus / jacobian * elastic_modulus;

if(no_damage) return SUANPAN_SUCCESS;
const auto dgamma = elastic_modulus / jacobian;
trial_stiffness -= dgamma * elastic_modulus;

const auto damage = compute_compression_damage(kc);
if(apply_damage) {
const auto damage = compute_compression_damage(kc);
const auto damage_factor = 1. - damage(0);

trial_stiffness *= 1. - damage(0);
trial_stiffness -= trial_stress * damage(1) * elastic_modulus / jacobian;
trial_stiffness *= damage_factor;
trial_stiffness -= trial_stress * damage(1) * dgamma;

trial_stress *= 1. - damage(0);
trial_stress *= damage_factor;
}

return SUANPAN_SUCCESS;
}
Expand All @@ -129,12 +141,13 @@ int NonlinearK4::compute_crack_close_branch() {

// account for entering
const auto net_strain = fabs(incre_strain(0)) - std::max(0., current_stress(0)) / elastic_modulus;
auto incre = net_strain * elastic_modulus / jacobian;
const auto dgamma = elastic_modulus / jacobian;
auto incre = net_strain * dgamma;

// physically, the tension plastic strain is the crack opening, closing the crack should not exceed the opening
// ensure the crack plastic strain is bounded by the tension plastic strain
if(incre > kt - kk) incre = kt - kk;
else trial_stiffness -= elastic_modulus / jacobian * elastic_modulus; // otherwise, the stiffness is degraded during the closing phase
else trial_stiffness -= dgamma * elastic_modulus; // otherwise, the stiffness is degraded during the closing phase

const auto incre_ep = incre * suanpan::sign(trial_stress(0));

Expand Down Expand Up @@ -173,7 +186,7 @@ int NonlinearK4::update_trial_status(const vec& n_strain) {
const auto& current_kt = current_history(1);
const auto& current_kk = current_history(3);

trial_stress = elastic_modulus * (trial_strain - plastic_strain);
trial_stress = (trial_stiffness = elastic_modulus) * (trial_strain - plastic_strain);

if(trial_stress(0) < 0. && incre_strain(0) < 0. && current_kt > current_kk && SUANPAN_SUCCESS != compute_crack_close_branch()) return SUANPAN_FAIL;

Expand Down
2 changes: 1 addition & 1 deletion Material/Material1D/Concrete/NonlinearK4.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -38,7 +38,7 @@ struct DataNonlinearK4 {
class NonlinearK4 : protected DataNonlinearK4, public Material1D {
static constexpr unsigned max_iteration = 20;

const bool no_damage = true;
const bool apply_damage = true;

[[nodiscard]] virtual vec2 compute_tension_backbone(double) const = 0;
[[nodiscard]] virtual vec2 compute_compression_backbone(double) const = 0;
Expand Down

0 comments on commit 2c1b789

Please sign in to comment.