compareAbaqusAnsysRotorEigenfrequencies.py

You can view and download this file on Github: compareAbaqusAnsysRotorEigenfrequencies.py

  1#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  2# This is an EXUDYN example
  3#
  4# Details:  Example to compute eigenfrequencies of a rotor
  5#           NOTE that this example requires files from a subdirectory testData as provided on github
  6#
  7# Author:   Stefan Holzinver, Johannes Gerstmayr
  8# Date:     2020-05-18
  9#
 10# Copyright:This file is part of Exudyn. Exudyn is free software. You can redistribute it and/or modify it under the terms of the Exudyn license. See 'LICENSE.txt' for more details.
 11#
 12#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 13
 14import numpy as np
 15from scipy.sparse import linalg
 16import scipy as sp
 17
 18import exudyn as exu
 19from exudyn.utilities import * #includes itemInterface and rigidBodyUtilities
 20import exudyn.graphics as graphics #only import if it does not conflict
 21from exudyn.FEM import *
 22
 23useGraphics = True #without test
 24#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 25#you can erase the following lines and all exudynTestGlobals related operations if this is not intended to be used as TestModel:
 26try: #only if called from test suite
 27    from modelUnitTests import exudynTestGlobals #for globally storing test results
 28    useGraphics = exudynTestGlobals.useGraphics
 29except:
 30    class ExudynTestGlobals:
 31        pass
 32    exudynTestGlobals = ExudynTestGlobals()
 33#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 34
 35numberOfModes = 18
 36useSparseSolverRoutine = False
 37
 38errorResult = 0
 39
 40testDataDir = "testData/"
 41#if useGraphics:
 42#    testDataDir = "testData/"
 43
 44###############################################################################
 45# Ansys - lumped mass matrix formulation - Sparse Matrix - MMF format
 46###############################################################################
 47
 48#read finite element model
 49exudynTestGlobals.testResult = 0
 50for testNumber in range(2):
 51    if testNumber == 1:
 52        useSparseSolverRoutine = True
 53
 54    fem = FEMinterface()
 55    #exu.Print("read matrices ...")
 56    fem.ReadMassMatrixFromAnsys(fileName=testDataDir + 'rotorAnsysMassMatrixSparse.txt',
 57                                dofMappingVectorFile=testDataDir + 'rotorAnsysMatrixDofMappingVector.txt')
 58    fem.ReadStiffnessMatrixFromAnsys(fileName=testDataDir + 'rotorAnsysStiffnessMatrixSparse.txt',
 59                                dofMappingVectorFile=testDataDir + 'rotorAnsysMatrixDofMappingVector.txt')
 60
 61    #exu.Print("compute eigenmodes ...")
 62    fem.ComputeEigenmodes(numberOfModes, useSparseSolver = useSparseSolverRoutine)
 63
 64    if useGraphics:
 65        exu.Print('natural frequencies from Ansys model (Lumped Mass Matrix, MMF-Format)', fem.GetEigenFrequenciesHz()[0:numberOfModes])
 66
 67    if not useSparseSolverRoutine:
 68        f6 = fem.GetEigenFrequenciesHz()[6] - 104.63701055079783 #reference solution
 69    else:
 70        #compare with dense matrix solution;
 71        #due to random initialization, the results change with every computation ==> approx. 1e-11 repeatability
 72        f6 = fem.GetEigenFrequenciesHz()[6] - 104.6370105507865
 73    f6*=1e-6 #use offset also for abaqus, as it gives non-reproducible results in dense case (32/64bit?)
 74    exu.Print('natural frequencies from Ansys model, sparse=',str(useSparseSolverRoutine),":", fem.GetEigenFrequenciesHz()[6] )
 75    errorResult += f6
 76
 77    exudynTestGlobals.testResult += 1e-6*fem.GetEigenFrequenciesHz()[6]
 78
 79    ###############################################################################
 80    # Abaqus
 81    ###############################################################################
 82
 83    #read finite element model
 84    fem = FEMinterface()
 85    #exu.Print("read matrices ...")
 86    fem.ReadMassMatrixFromAbaqus(fileName=testDataDir + 'rotorDiscTestMASS1.mtx')
 87    fem.ReadStiffnessMatrixFromAbaqus(fileName=testDataDir + 'rotorDiscTestSTIF1.mtx')
 88    fem.ScaleStiffnessMatrix(1e-2) #in Ansys, the stiffness matrix is already scaled!
 89
 90    #exu.Print("compute eigenmodes ...")
 91    fem.ComputeEigenmodes(numberOfModes, useSparseSolver = useSparseSolverRoutine)
 92
 93    if useGraphics:
 94        exu.Print('natural frequencies from Abaqus model (Lumped Mass Matrix)',fem.GetEigenFrequenciesHz()[0:numberOfModes])
 95
 96    if not useSparseSolverRoutine:
 97        f6 = fem.GetEigenFrequenciesHz()[6] - 104.6370132606291 #reference solution
 98    else:
 99        #compare with dense matrix solution;
100        #due to random initialization, the results change with every computation ==> approx. 1e-11 repeatability
101        f6 = fem.GetEigenFrequenciesHz()[6] - 104.6370132606291
102    f6*=1e-6 #use offset also for abaqus, as it gives non-reproducible results in dense case (32/64bit?)
103    exu.Print('natural frequencies from Abaqus model, sparse=',str(useSparseSolverRoutine),":", fem.GetEigenFrequenciesHz()[6] )
104    errorResult += f6
105    exudynTestGlobals.testResult += 1e-6*fem.GetEigenFrequenciesHz()[6]
106
107exu.Print('error of compareAbaqusAnsysRotorEigenfrequencies (due to sparse solver)=',errorResult)
108if abs(errorResult) < 1e-15: #usually of size 1e-17
109    errorResult = 0 #due to randomized sparse solver results, take this treshold!
110
111exu.Print('solution of compareAbaqusAnsysRotorEigenfrequencies (with treshold)=',errorResult)
112exudynTestGlobals.testError = errorResult #2020-05-22: 0
113#exudynTestGlobals.testResult computed above