Shell strain principal and von Mises difference with MSC

[victorkemp]

The strain components are the same, but Mystran has 9.01388E-04 and 5.00000E-04 for von Mises where MSC has 3.333333E-04. Principal strains are also different.

Affects v18.0.0 and current dev branch.

Input

SOL 101
CEND
DISPLACEMENT = ALL
STRAIN(FIBER, CENTER) = ALL
SPC = 1
LOAD = 1
BEGIN BULK
GRID           1       0   0.0+0   0.0+0   0.0+0       0
GRID           2       0   1.5+0   0.0+0   0.0+0       0
GRID           3       0   0.0+0   7.5-1   0.0+0       0
GRID           4       0   7.5-1   7.5-1   0.0+0       0
GRID           5       0   7.5-1   0.0+0   0.0+0       0
GRID           6       0   1.5+0   7.5-1   0.0+0       0
GRID           7       0   0.0+0   1.0+0   0.0+0       0
GRID           8       0   1.5+0   1.0+0   0.0+0       0
GRID           9       0   0.0+0  1.75+0   0.0+0       0
GRID          10       0   7.5-1  1.75+0   0.0+0       0
GRID          11       0   7.5-1   1.0+0   0.0+0       0
GRID          12       0   1.5+0  1.75+0   0.0+0       0
MAT1           1  2.0+11  1.0+11   0.0+0   0.0+0   0.0+0   0.0+0
SPC            1       1       2   0.0+0
SPC            1       1       1   0.0+0
SPC            1       1       3   0.0+0
SPC            1       1       5   0.0+0
SPC            1       3       1   0.0+0
SPC            1       3       3   0.0+0
SPC            1       3       5   0.0+0
SPC            1       7       2   0.0+0
SPC            1       7       1   0.0+0
SPC            1       7       3   0.0+0
SPC            1       7       5   0.0+0
SPC            1       9       1   0.0+0
SPC            1       9       3   0.0+0
SPC            1       9       5   0.0+0
PSHELL         1       1   4.0-2       1               1           0.0+0
CQUAD4         1       1       5       2       6       4   0.0+0   0.0+0
CQUAD4         2       1       1       5       4       3   0.0+0   0.0+0
PSHELL         2       1   4.0-2       1               1           0.0+0
CTRIA3         3       2      11      12      10   0.0+0   0.0+0
CTRIA3         4       2       7      10       9   0.0+0   0.0+0
CTRIA3         5       2      11       8      12   0.0+0   0.0+0
CTRIA3         6       2       7      11      10   0.0+0   0.0+0
MOMENT         2       6       0   1.0+0   0.0+0  -1.0+4   0.0+0
MOMENT         2       2       0   1.0+0   0.0+0  -1.0+4   0.0+0
MOMENT         2      12       0   1.0+0   0.0+0  -1.0+4   0.0+0
MOMENT         2       8       0   1.0+0   0.0+0  -1.0+4   0.0+0
LOAD           1   1.0+0   1.0+0       2
ENDDATA 

MSC f06:

                             S T R A I N S   I N   T R I A N G U L A R   E L E M E N T S   ( T R I A 3 )
  ELEMENT      FIBER                STRAINS IN ELEMENT COORD SYSTEM             PRINCIPAL  STRAINS (ZERO SHEAR)                 
    ID.       DISTANCE           NORMAL-X       NORMAL-Y      SHEAR-XY       ANGLE         MAJOR           MINOR        VON MISES
0       3   -2.000000E-02      2.500000E-04   2.500000E-04  -5.000000E-04   -45.0000    5.000000E-04   -7.047314E-19   3.333333E-04
             2.000000E-02     -2.500000E-04  -2.500000E-04   5.000000E-04    45.0000    7.047314E-19   -5.000000E-04   3.333333E-04
0       4   -2.000000E-02      2.500000E-04   2.500000E-04  -5.000000E-04   -45.0000    5.000000E-04    7.047314E-19   3.333333E-04
             2.000000E-02     -2.500000E-04  -2.500000E-04   5.000000E-04    45.0000   -7.047314E-19   -5.000000E-04   3.333333E-04
0       5   -2.000000E-02      5.000000E-04  -2.405483E-18  -2.775558E-19    -0.0000    5.000000E-04   -2.385245E-18   3.333333E-04
             2.000000E-02     -5.000000E-04   2.405483E-18   2.775558E-19    90.0000    2.385245E-18   -5.000000E-04   3.333333E-04
0       6   -2.000000E-02      5.000000E-04  -6.649773E-19  -1.249001E-18    -0.0000    5.000000E-04   -6.505213E-19   3.333333E-04
             2.000000E-02     -5.000000E-04   6.649773E-19   1.249001E-18    90.0000    6.505213E-19   -5.000000E-04   3.333333E-04

Mystran f06

                                     E L E M E N T   S T R A I N S   I N   L O C A L   E L E M E N T   C O O R D I N A T E   S Y S T E M
                                                                 F O R   E L E M E N T   T Y P E   T R I A 3  
  Element    Location      Fibre        Strains In Element Coord System       Principal Strains (Zero Shear)                 Transverse   Transverse
     ID                   Distance     Normal-X     Normal-Y      Shear-XY     Angle     Major        Minor      von Mises    Shear-XZ     Shear-YZ 

        3    Anywhere   -2.00000E-02  2.50000E-04  2.50000E-04 -5.00000E-04  -45.000  7.50000E-04 -2.50000E-04  9.01388E-04  3.90313E-18  8.67362E-19
             in elem     2.00000E-02 -2.50000E-04 -2.50000E-04  5.00000E-04   45.000  2.50000E-04 -7.50000E-04  9.01388E-04

        4    Anywhere   -2.00000E-02  2.50000E-04  2.50000E-04 -5.00000E-04  -45.000  7.50000E-04 -2.50000E-04  9.01388E-04  1.92151E-18 -6.37855E-18
             in elem     2.00000E-02 -2.50000E-04 -2.50000E-04  5.00000E-04   45.000  2.50000E-04 -7.50000E-04  9.01388E-04

        5    Anywhere   -2.00000E-02  5.00000E-04  6.17897E-19 -1.38778E-18   -0.000  5.00000E-04  6.50521E-19  5.00000E-04 -6.07153E-18  1.73472E-18
             in elem     2.00000E-02 -5.00000E-04 -6.17897E-19  1.38778E-18   90.000 -6.50521E-19 -5.00000E-04  5.00000E-04

        6    Anywhere   -2.00000E-02  5.00000E-04 -2.29227E-20  1.31839E-18    0.000  5.00000E-04  0.00000E+00  5.00000E-04 -6.07153E-18 -4.33681E-19
             in elem     2.00000E-02 -5.00000E-04  2.29227E-20 -1.31839E-18  -90.000  0.00000E+00 -5.00000E-04  5.00000E-04

[SteveDoyle2]

The equations for von mises and principal stresses/strains are:

def von_mises_2d(oxx: np.ndarray,
                 oyy: np.ndarray,
                 txy: np.ndarray,
                 is_stress: bool) -> np.ndarray:
    """verified for stress/strain"""
    if is_stress:
        ovm = np.sqrt(oxx**2 + oyy**2 - oxx*oyy + 3*(txy**2) )
    else:
        # center = (oxx + oyy) / 2
        # radius = np.sqrt((oxx - oyy) ** 2 / 4 + txy ** 2 / 4)
        # eig_max = center + radius
        # eig_min = center - radius
        # return von_mises_2d(eig_max, eig_min, txy*0, is_stress=True)
        ovm = np.sqrt(4/9 * (oxx**2 + oyy**2 - oxx*oyy) + (txy**2)/3 )
    return ovm

def principal_2d(oxx: np.ndarray,
                 oyy: np.ndarray,
                 txy: np.ndarray,
                 is_stress: bool) -> tuple[np.ndarray, np.ndarray]:
    """
    https://www.simscale.com/docs/simwiki/fea-finite-element-analysis/principal-stress-and-principal-strain/
    """
    if is_stress:
        center = (oxx + oyy) / 2
        radius = np.sqrt((oxx - oyy) ** 2 / 4 + txy ** 2)
    else:
        center = (oxx + oyy) / 2
        radius = 0.5 * np.sqrt((oxx - oyy) ** 2 + txy ** 2)
    eig_max = center + radius
    eig_min = center - radius
    return eig_max, eig_min