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 *
 20from exudyn.FEM import *
 21
 22useGraphics = True #without test
 23#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 24#you can erase the following lines and all exudynTestGlobals related operations if this is not intended to be used as TestModel:
 25try: #only if called from test suite
 26    from modelUnitTests import exudynTestGlobals #for globally storing test results
 27    useGraphics = exudynTestGlobals.useGraphics
 28except:
 29    class ExudynTestGlobals:
 30        pass
 31    exudynTestGlobals = ExudynTestGlobals()
 32#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 33
 34numberOfModes = 18
 35useSparseSolverRoutine = False
 36
 37errorResult = 0
 38
 39testDataDir = "testData/"
 40#if useGraphics:
 41#    testDataDir = "testData/"
 42
 43###############################################################################
 44# Ansys - lumped mass matrix formulation - Sparse Matrix - MMF format
 45###############################################################################
 46
 47#read finite element model
 48exudynTestGlobals.testResult = 0
 49for testNumber in range(2):
 50    if testNumber == 1:
 51        useSparseSolverRoutine = True
 52
 53    fem = FEMinterface()
 54    #exu.Print("read matrices ...")
 55    fem.ReadMassMatrixFromAnsys(fileName=testDataDir + 'rotorAnsysMassMatrixSparse.txt',
 56                                dofMappingVectorFile=testDataDir + 'rotorAnsysMatrixDofMappingVector.txt')
 57    fem.ReadStiffnessMatrixFromAnsys(fileName=testDataDir + 'rotorAnsysStiffnessMatrixSparse.txt',
 58                                dofMappingVectorFile=testDataDir + 'rotorAnsysMatrixDofMappingVector.txt')
 59
 60    #exu.Print("compute eigenmodes ...")
 61    fem.ComputeEigenmodes(numberOfModes, useSparseSolver = useSparseSolverRoutine)
 62
 63    if useGraphics:
 64        exu.Print('natural frequencies from Ansys model (Lumped Mass Matrix, MMF-Format)', fem.GetEigenFrequenciesHz()[0:numberOfModes])
 65
 66    if not useSparseSolverRoutine:
 67        f6 = fem.GetEigenFrequenciesHz()[6] - 104.63701055079783 #reference solution
 68    else:
 69        #compare with dense matrix solution;
 70        #due to random initialization, the results change with every computation ==> approx. 1e-11 repeatability
 71        f6 = fem.GetEigenFrequenciesHz()[6] - 104.6370105507865
 72    f6*=1e-6 #use offset also for abaqus, as it gives non-reproducible results in dense case (32/64bit?)
 73    exu.Print('natural frequencies from Ansys model, sparse=',str(useSparseSolverRoutine),":", fem.GetEigenFrequenciesHz()[6] )
 74    errorResult += f6
 75
 76    exudynTestGlobals.testResult += 1e-6*fem.GetEigenFrequenciesHz()[6]
 77
 78    ###############################################################################
 79    # Abaqus
 80    ###############################################################################
 81
 82    #read finite element model
 83    fem = FEMinterface()
 84    #exu.Print("read matrices ...")
 85    fem.ReadMassMatrixFromAbaqus(fileName=testDataDir + 'rotorDiscTestMASS1.mtx')
 86    fem.ReadStiffnessMatrixFromAbaqus(fileName=testDataDir + 'rotorDiscTestSTIF1.mtx')
 87    fem.ScaleStiffnessMatrix(1e-2) #in Ansys, the stiffness matrix is already scaled!
 88
 89    #exu.Print("compute eigenmodes ...")
 90    fem.ComputeEigenmodes(numberOfModes, useSparseSolver = useSparseSolverRoutine)
 91
 92    if useGraphics:
 93        exu.Print('natural frequencies from Abaqus model (Lumped Mass Matrix)',fem.GetEigenFrequenciesHz()[0:numberOfModes])
 94
 95    if not useSparseSolverRoutine:
 96        f6 = fem.GetEigenFrequenciesHz()[6] - 104.6370132606291 #reference solution
 97    else:
 98        #compare with dense matrix solution;
 99        #due to random initialization, the results change with every computation ==> approx. 1e-11 repeatability
100        f6 = fem.GetEigenFrequenciesHz()[6] - 104.6370132606291
101    f6*=1e-6 #use offset also for abaqus, as it gives non-reproducible results in dense case (32/64bit?)
102    exu.Print('natural frequencies from Abaqus model, sparse=',str(useSparseSolverRoutine),":", fem.GetEigenFrequenciesHz()[6] )
103    errorResult += f6
104    exudynTestGlobals.testResult += 1e-6*fem.GetEigenFrequenciesHz()[6]
105
106exu.Print('error of compareAbaqusAnsysRotorEigenfrequencies (due to sparse solver)=',errorResult)
107if abs(errorResult) < 1e-15: #usually of size 1e-17
108    errorResult = 0 #due to randomized sparse solver results, take this treshold!
109
110exu.Print('solution of compareAbaqusAnsysRotorEigenfrequencies (with treshold)=',errorResult)
111exudynTestGlobals.testError = errorResult #2020-05-22: 0
112#exudynTestGlobals.testResult computed above