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