diff --git a/src/constitutive/TACSBladeStiffenedShellConstitutive.cpp b/src/constitutive/TACSBladeStiffenedShellConstitutive.cpp index bc6c02630..7825609e8 100644 --- a/src/constitutive/TACSBladeStiffenedShellConstitutive.cpp +++ b/src/constitutive/TACSBladeStiffenedShellConstitutive.cpp @@ -1225,6 +1225,20 @@ TacsScalar TACSBladeStiffenedShellConstitutive::evalDesignFieldValue( return 0.0; } +TacsScalar TACSBladeStiffenedShellConstitutive::evalFailureFieldValue( + int elemIndex, const double pt[], const TacsScalar X[], + const TacsScalar strain[], int failIndex) { + if (failIndex == 0) { + return this->evalFailure(elemIndex, pt, X, strain); + } else if (failIndex >= 1 && failIndex <= this->NUM_FAILURES) { + TacsScalar fails[this->NUM_FAILURES]; + computeFailureValues(strain, fails); + return fails[failIndex - 1]; + } else { + return 0.0; + } +} + TacsScalar TACSBladeStiffenedShellConstitutive::computeEffectiveBendingThickness() { TacsScalar IStiff = this->computeStiffenerIzz(); diff --git a/src/constitutive/TACSBladeStiffenedShellConstitutive.h b/src/constitutive/TACSBladeStiffenedShellConstitutive.h index 0e1d93148..4e466edc5 100644 --- a/src/constitutive/TACSBladeStiffenedShellConstitutive.h +++ b/src/constitutive/TACSBladeStiffenedShellConstitutive.h @@ -318,6 +318,11 @@ class TACSBladeStiffenedShellConstitutive : public TACSShellConstitutive { TacsScalar evalDesignFieldValue(int elemIndex, const double pt[], const TacsScalar X[], int index); + // Retrieve failure mode values for plotting purposes + TacsScalar evalFailureFieldValue(int elemIndex, const double pt[], + const TacsScalar X[], + const TacsScalar strain[], int failIndex); + /** * @brief Compute the effective tensile thickness of the stiffened shell * diff --git a/src/constitutive/TACSConstitutive.h b/src/constitutive/TACSConstitutive.h index 80107f98d..d04545aae 100644 --- a/src/constitutive/TACSConstitutive.h +++ b/src/constitutive/TACSConstitutive.h @@ -578,6 +578,26 @@ class TACSConstitutive : public TACSObject { return 0.0; } + /** + Evaluate the failure index (if defined) at the given quadrature point + + @param elemIndex The local element index + @param pt The parametric point + @param X The physical node location + @param strain The strain value + @param The failure index + @return The failure value + */ + virtual TacsScalar evalFailureFieldValue(int elemIndex, const double pt[], + const TacsScalar X[], + const TacsScalar strain[], + int failIndex) { + if (failIndex == 0) { + return evalFailure(elemIndex, pt, X, strain); + } + return 0.0; + } + /** Compute a two-dimensional representation of the failure envelope */ diff --git a/src/elements/TACSElementTypes.cpp b/src/elements/TACSElementTypes.cpp index 6269b5d1c..76b4293b6 100644 --- a/src/elements/TACSElementTypes.cpp +++ b/src/elements/TACSElementTypes.cpp @@ -85,7 +85,7 @@ int TacsGetOutputComponentCount(ElementType etype, int comp) { } else if (comp == TACS_OUTPUT_STRESSES) { return 9; } else if (comp == TACS_OUTPUT_EXTRAS) { - return 8; + return 13; } else if (comp == TACS_OUTPUT_LOADS) { return 6; } @@ -298,20 +298,30 @@ const char *TacsGetOutputComponentName(ElementType etype, int comp, int index) { } else if (comp == TACS_OUTPUT_EXTRAS) { switch (index) { case 0: - return "failure"; + return "failure0"; case 1: - return "dv1"; + return "failure1"; case 2: - return "dv2"; + return "failure2"; case 3: - return "dv3"; + return "failure3"; case 4: - return "dv4"; + return "failure4"; case 5: - return "dv5"; + return "failure5"; case 6: - return "dv6"; + return "dv1"; case 7: + return "dv2"; + case 8: + return "dv3"; + case 9: + return "dv4"; + case 10: + return "dv5"; + case 11: + return "dv6"; + case 12: return "dv7"; default: return NULL; diff --git a/src/elements/shell/TACSBeamElement.h b/src/elements/shell/TACSBeamElement.h index 26cd69a6e..35ab9f76f 100644 --- a/src/elements/shell/TACSBeamElement.h +++ b/src/elements/shell/TACSBeamElement.h @@ -2222,11 +2222,16 @@ void TACSBeamElement::getOutputData( data += 9; } if (write_flag & TACS_OUTPUT_EXTRAS) { - data[0] = con->evalFailure(elemIndex, pt, X0.x, e); - data[1] = con->evalDesignFieldValue(elemIndex, pt, X0.x, 0); - data[2] = con->evalDesignFieldValue(elemIndex, pt, X0.x, 1); - data[3] = con->evalDesignFieldValue(elemIndex, pt, X0.x, 2); - data += 4; + data[0] = con->evalFailureFieldValue(elemIndex, pt, X0.x, e, 0); + data[1] = con->evalFailureFieldValue(elemIndex, pt, X0.x, e, 1); + data[2] = con->evalFailureFieldValue(elemIndex, pt, X0.x, e, 2); + data[3] = con->evalFailureFieldValue(elemIndex, pt, X0.x, e, 3); + data[4] = con->evalFailureFieldValue(elemIndex, pt, X0.x, e, 4); + data[5] = con->evalFailureFieldValue(elemIndex, pt, X0.x, e, 5); + data[6] = con->evalDesignFieldValue(elemIndex, pt, X0.x, 0); + data[7] = con->evalDesignFieldValue(elemIndex, pt, X0.x, 1); + data[8] = con->evalDesignFieldValue(elemIndex, pt, X0.x, 2); + data += 9; } } } diff --git a/src/elements/shell/TACSShellElement.h b/src/elements/shell/TACSShellElement.h index 6bf3ddf77..06e3dc090 100644 --- a/src/elements/shell/TACSShellElement.h +++ b/src/elements/shell/TACSShellElement.h @@ -1487,15 +1487,20 @@ void TACSShellElement::getOutputData( data += 9; } if (write_flag & TACS_OUTPUT_EXTRAS) { - data[0] = con->evalFailure(elemIndex, pt, X, e); - data[1] = con->evalDesignFieldValue(elemIndex, pt, X, 0); - data[2] = con->evalDesignFieldValue(elemIndex, pt, X, 1); - data[3] = con->evalDesignFieldValue(elemIndex, pt, X, 2); - data[4] = con->evalDesignFieldValue(elemIndex, pt, X, 3); - data[5] = con->evalDesignFieldValue(elemIndex, pt, X, 4); - data[6] = con->evalDesignFieldValue(elemIndex, pt, X, 5); - data[7] = con->evalDesignFieldValue(elemIndex, pt, X, 6); - data += 8; + data[0] = con->evalFailureFieldValue(elemIndex, pt, X, e, 0); + data[1] = con->evalFailureFieldValue(elemIndex, pt, X, e, 1); + data[2] = con->evalFailureFieldValue(elemIndex, pt, X, e, 2); + data[3] = con->evalFailureFieldValue(elemIndex, pt, X, e, 3); + data[4] = con->evalFailureFieldValue(elemIndex, pt, X, e, 4); + data[5] = con->evalFailureFieldValue(elemIndex, pt, X, e, 5); + data[6] = con->evalDesignFieldValue(elemIndex, pt, X, 0); + data[7] = con->evalDesignFieldValue(elemIndex, pt, X, 1); + data[8] = con->evalDesignFieldValue(elemIndex, pt, X, 2); + data[9] = con->evalDesignFieldValue(elemIndex, pt, X, 3); + data[10] = con->evalDesignFieldValue(elemIndex, pt, X, 4); + data[11] = con->evalDesignFieldValue(elemIndex, pt, X, 5); + data[12] = con->evalDesignFieldValue(elemIndex, pt, X, 6); + data += 13; } } } diff --git a/src/elements/shell/TACSThermalShellElement.h b/src/elements/shell/TACSThermalShellElement.h index 12406156b..049ce7017 100644 --- a/src/elements/shell/TACSThermalShellElement.h +++ b/src/elements/shell/TACSThermalShellElement.h @@ -1451,11 +1451,16 @@ void TACSThermalShellElement::getOutputData( data += 9; } if (write_flag & TACS_OUTPUT_EXTRAS) { - data[0] = con->evalFailure(elemIndex, pt, X, e); - data[1] = con->evalDesignFieldValue(elemIndex, pt, X, 0); - data[2] = con->evalDesignFieldValue(elemIndex, pt, X, 1); - data[3] = con->evalDesignFieldValue(elemIndex, pt, X, 2); - data += 4; + data[0] = con->evalFailureFieldValue(elemIndex, pt, X, e, 0); + data[1] = con->evalFailureFieldValue(elemIndex, pt, X, e, 1); + data[2] = con->evalFailureFieldValue(elemIndex, pt, X, e, 2); + data[3] = con->evalFailureFieldValue(elemIndex, pt, X, e, 3); + data[4] = con->evalFailureFieldValue(elemIndex, pt, X, e, 4); + data[5] = con->evalFailureFieldValue(elemIndex, pt, X, e, 5); + data[6] = con->evalDesignFieldValue(elemIndex, pt, X, 0); + data[7] = con->evalDesignFieldValue(elemIndex, pt, X, 1); + data[8] = con->evalDesignFieldValue(elemIndex, pt, X, 2); + data += 9; } } }