5
5
Harmonic analysis using the frequency-sweep Krylov method
6
6
=========================================================
7
7
8
- This example shows how to use the Krylov method in `PyMAPDL <https://mapdl.docs.pyansys.com/ >`_ for harmonic
9
- analysis.
8
+ This example shows how to use the frequency-sweep Krylov method
9
+ implemented in PyMAPDL. For more information, including the
10
+ theory behind this method, see `Frequency-Sweep Harmonic Analysis via the Krylov Method
11
+ <ansys_krylov_sweep_harmonic_analysis_> `_
12
+ in the **Structural Analysis ** guide for Mechanical APDL.
10
13
11
- These are the main steps required:
14
+ Overview
15
+ --------
12
16
13
- - Generate a Krylov subspace for model reduction in the harmonic analysis using
14
- the :func: ` KrylovSolver.gensubspace() <ansys.mapdl.core.krylov.KrylovSolver.gensubspace> `
15
- method .
17
+ This example uses the frequency-sweep Krylov method to perform a harmonic analysis
18
+ on a cylindrical acoustic duct and study the response of the system over
19
+ a range of frequencies .
16
20
17
- - Reduce the system of equations and solve at each frequency using the
18
- :func: ` KrylovSolver.solve() <ansys.mapdl.core.krylov.KrylovSolver.solve> ` method .
21
+ The model is a cylindrical acoustic duct with pressure load on one end
22
+ and output impedance on the other end .
19
23
20
- - Expand the reduced solution back to the FE space using the
21
- :func: `KrylovSolver.expand() <ansys.mapdl.core.krylov.KrylovSolver.expand> ` method.
24
+ These are the main steps required:
22
25
23
- Problem Description
24
- -------------------
26
+ - Use the :func: ` KrylovSolver.gensubspace() <ansys.mapdl.core.krylov.KrylovSolver.gensubspace> `
27
+ method to generate a Krylov subspace for model reduction in the harmonic analysis.
25
28
26
- To perform an harmonic analysis on a cylindrical acoustic duct using the Krylov method
27
- and study the response of the system over a range of frequencies .
29
+ - Use the :func: ` KrylovSolver.solve() <ansys.mapdl.core.krylov.KrylovSolver.solve> `
30
+ method to reduce the system of equations and solve at each frequency .
28
31
32
+ - Use the :func: `KrylovSolver.expand() <ansys.mapdl.core.krylov.KrylovSolver.expand> ` method
33
+ to expand the reduced solution back to the FE space.
29
34
30
- The model is a cylindrical acoustic duct with pressure load on one end
31
- and output impedance on the other end.
35
+ Perform required imports
36
+ ------------------------
37
+ Perform required imports and launch mapdl.
32
38
33
39
.. code :: ipython3
34
40
@@ -42,16 +48,14 @@ and output impedance on the other end.
42
48
mapdl.clear()
43
49
mm = mapdl.math
44
50
51
+
52
+ Define parameters
53
+ -----------------
45
54
46
-
47
-
48
- Parameters definition
49
- ~~~~~~~~~~~~~~~~~~~~~
50
-
51
- This section defines some geometry parameters and analysis settings.
52
- As mentioned earlier, the geometry is a cylinder defined by its radius (``cyl_r ``) and its length (``cyl_L ``).
53
- The length of the duct is such that three complete wavelengths (``no_wl ``) can fit in its length and can have
54
- ten elements per wavelength.
55
+ Define some geometry parameters and analysis settings. As mentioned earlier, the geometry
56
+ is a cylinder defined by its radius (``cyl_r ``) and its length (``cyl_L ``). The length
57
+ of the duct is such that three complete wavelengths (``no_wl ``) can fit in its length
58
+ and can have ten elements per wavelength.
55
59
56
60
.. code :: ipython3
57
61
@@ -75,11 +79,10 @@ ten elements per wavelength.
75
79
nelem_wl = 10 # no of elements per wavelength
76
80
tol_elem = nelem_wl * no_wl # total number of elements across length
77
81
78
- Element and material definition
79
- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
80
-
81
- Here we will assign fluid medium (Air) properties to the duct.
82
- We will use Fluid 220 (Keyopt(2)=1) with one degree of freedom per node : pressure,
82
+ Define element type and materials
83
+ ---------------------------------
84
+ Assign fluid medium (air) properties to the duct. This example
85
+ uses Fluid 220 (``Keyopt(2)=1 ``) with one degree of freedom per node (pressure),
83
86
with no FSI interface in the element.
84
87
85
88
.. code :: ipython3
@@ -92,44 +95,45 @@ with no FSI interface in the element.
92
95
mapdl.mp("VISC", 1, visco)
93
96
94
97
95
- Geometry definition
96
- ~~~~~~~~~~~~~~~~~~~
98
+ Define geometry
99
+ ---------------
97
100
98
- Here we discuss creating a cylinder of required dimensions and splitting it into
99
- 4 segments for uniform generation of the mesh in each segment.
101
+ Create a cylinder of the required dimensions and split it into
102
+ four segments for uniform generation of the mesh in each segment.
100
103
101
104
.. code :: ipython3
102
105
103
- # Setting back to default
106
+ # Set back to default
104
107
mapdl.csys(0)
105
108
106
- # Rotating working plane for the cylinder generation
109
+ # Rotate working plane for the cylinder generation
107
110
mapdl.wpcsys(-1)
108
111
mapdl.wprota(thzx=90)
109
112
110
- # Generating a circular area with specified radius
113
+ # Generate a circular area with a specified radius
111
114
mapdl.cyl4(0, 0, cyl_r)
112
115
113
116
mapdl.wpcsys(-1)
114
117
115
- # Extrude the circular area to generate cylinder of specified length
118
+ # Extrude the circular area to generate a cylinder of specified length
116
119
mapdl.vext("ALL", dx=cyl_L)
117
120
118
- # Segment the cylinder into four quadrants to create a more uniform mesh
121
+ # Split the cylinder into four segments to create a more uniform mesh
119
122
mapdl.vsbw("ALL", keep='DELETE')
120
123
mapdl.wprota(thzx=90)
121
- mapdl.vsbw("ALL", keep='DELETE') # Why this is needed?
124
+ mapdl.vsbw("ALL", keep='DELETE')
122
125
123
126
mapdl.wpcsys(-1)
124
127
125
- # Creating a component with the created volume
128
+ # Create a component with the created volume
126
129
mapdl.cm('cm1', 'volu')
127
130
128
131
129
132
130
- Create mesh:
133
+ Create mesh
134
+ -----------
131
135
132
- Mesh the
136
+ Create the mesh and plot the FE model.
133
137
134
138
.. code :: ipython3
135
139
@@ -154,28 +158,23 @@ a length element size constraint to the longitudinal lines.
154
158
mapdl.lsel("U", 'loc', 'x', 0)
155
159
mapdl.lsel("U", 'loc', 'x', cyl_L)
156
160
157
- # Apply length constraint.
161
+ # Apply length constraint
158
162
mapdl.lesize('ALL',ndiv = tol_elem)
159
163
mapdl.lsla()
160
164
161
165
# Mesh
162
166
mapdl.vsweep('ALL')
163
-
164
167
mapdl.allsel()
165
-
166
-
167
-
168
- Plot FE model:
169
-
170
- .. code :: ipython3
171
-
168
+
169
+ # Plot the FE model
172
170
mapdl.eplot()
173
171
172
+
174
173
.. image :: ../../../examples/extended_examples/Krylov/Harmonic_Analysis_using_krylov_pymapdl_files/Harmonic_Analysis_using_krylov_pymapdl_15_1.png
175
174
176
175
177
176
Define boundary conditions
178
- ~~~~~~~~~~~~~~~~~~~~~~~~~~
177
+ --------------------------
179
178
180
179
Apply pressure load on one end and output impedance on other end of the acoustic duct.
181
180
@@ -203,11 +202,10 @@ Apply pressure load on one end and output impedance on other end of the acoustic
203
202
mapdl.allsel()
204
203
205
204
206
-
207
205
Perform modal analysis
208
206
----------------------
209
207
210
- Get the first 10 natural frequency modes of the acoustic duct
208
+ Get the first 10 natural frequency modes of the acoustic duct.
211
209
212
210
.. code :: ipython3
213
211
@@ -243,28 +241,24 @@ Get the first 10 natural frequency modes of the acoustic duct
243
241
244
242
Run harmonic analysis using Krylov method
245
243
-----------------------------------------
244
+ Perform the following steps to run the harmonic analysis using the
245
+ frequency-sweep Krylov method.
246
246
247
- **Step 1 **: Generate full file.
247
+ **Step 1 **: Generate FULL file and initialize the `` Krylov `` class object .
248
248
249
249
.. code :: ipython3
250
250
251
251
mapdl.run('/SOLU')
252
- mapdl.antype('HARMIC') # HARMONIC ANALYSIS
252
+ mapdl.antype('HARMIC') # Set options for harmonic analysis
253
253
mapdl.hropt('KRYLOV')
254
254
mapdl.eqslv('SPARSE')
255
- mapdl.harfrq(0,1000) # set beginning and ending frequency
256
- mapdl.nsubst(100) # set the number of frequency increments
257
- mapdl.wrfull(1) # GENERATE . FULL FILE AND STOP
255
+ mapdl.harfrq(0,1000) # Set beginning and ending frequency
256
+ mapdl.nsubst(100) # Set the number of frequency increments
257
+ mapdl.wrfull(1) # Generate FULL file and stop
258
258
mapdl.solve()
259
259
mapdl.finish()
260
260
261
-
262
-
263
- Initialize Krylov class object.
264
-
265
- .. code :: ipython3
266
-
267
- dd = mapdl.krylov
261
+ dd = mapdl.krylov # Initialize Krylov class object
268
262
269
263
**Step 2 **: Generate a Krylov subspace of size/dimension 10 at frequency
270
264
500 Hz for model reduction.
@@ -273,7 +267,8 @@ Initialize Krylov class object.
273
267
274
268
Qz = dd.gensubspace(10, 500, check_orthogonality=True)
275
269
276
- The shape of the generated subspace.
270
+
271
+ Obtain the shape of the generated subspace.
277
272
278
273
.. code :: ipython3
279
274
@@ -292,7 +287,7 @@ from 0 Hz to 1000 Hz with ramped loading.
292
287
293
288
Yz = dd.solve(0, 1000, 100, ramped_load=True)
294
289
295
- The shape of the reduced solution generated is .
290
+ Obtain the shape of the reduced solution generated.
296
291
297
292
.. code :: ipython3
298
293
@@ -308,9 +303,9 @@ The shape of the reduced solution generated is.
308
303
309
304
.. code :: ipython3
310
305
311
- result = dd.expand(residual_computation=True, residual_algorithm="l2")
306
+ result = dd.expand(residual_computation=True, residual_algorithm="l2", return_solution = True )
312
307
313
- Results: Pressure distribution as a function of length
308
+ Plot the pressure distribution as a function of length
314
309
------------------------------------------------------
315
310
316
311
Plot the pressure distribution over the length of the duct on nodes where Y, Z coordinates are zero.
@@ -334,13 +329,13 @@ Load the last result substep to get the pressure for each of the selected nodes.
334
329
substep_index = 99
335
330
336
331
def get_pressure_at(node, step=1):
337
- """Get the pressure at a given node at a given step (by default first step)"""
332
+ """Get pressure at a given node at a given step (by default first step)"""
338
333
index_num = np.where(result[step]['node'] == node)
339
334
return result[step][index_num]
340
335
341
336
for each_node, loc in zip(ind, coords):
342
337
# Get pressure at the node
343
- pressure = get_pressure (each_node, substep_index)['x'][0]
338
+ pressure = get_pressure_at (each_node, substep_index)['x'][0]
344
339
345
340
#Calculate amplitude at 60 deg
346
341
magnitude = abs(pressure)
@@ -351,13 +346,13 @@ Load the last result substep to get the pressure for each of the selected nodes.
351
346
x_data.append(loc[0]) # X-Coordenate
352
347
y_data.append(pressure_a) # Nodal pressure at 60 degrees
353
348
354
- Sort the results according to the X- coordinate:
349
+ Sort the results according to the X coordinate.
355
350
356
351
.. code :: ipython3
357
352
358
353
sorted_x_data, sorted_y_data = zip(*sorted(zip(x_data, y_data)))
359
354
360
- Plot the calculated data:
355
+ Plot the calculated data.
361
356
362
357
.. code :: ipython3
363
358
@@ -378,11 +373,11 @@ Plot the calculated data:
378
373
.. image :: ../../../examples/extended_examples/Krylov/Harmonic_Analysis_using_krylov_pymapdl_files/Harmonic_Analysis_using_krylov_pymapdl_36_1.png
379
374
380
375
381
- Results: Plot frequency response function
382
- ------------------------------------------
376
+ Plot the frequency response function
377
+ ------------------------------------
383
378
384
379
Plot the frequency response function of any node along the length of the cylindrical duct.
385
- For example, let us plot the frequency response function for a node along 0.2 in X- direction of the duct.
380
+ This code plots the frequency response function for a node along 0.2 in the X direction of the duct.
386
381
387
382
.. code :: ipython3
388
383
@@ -391,7 +386,7 @@ For example, let us plot the frequency response function for a node along 0.2 in
391
386
392
387
393
388
Get the response of the system for the selected node
394
- over a range of frequencies [0- 1000 Hz]:
389
+ over a range of frequencies, such as 0 to 1000 Hz.
395
390
396
391
.. code :: python3
397
392
@@ -402,13 +397,13 @@ over a range of frequencies [0-1000 Hz]:
402
397
dic = {}
403
398
404
399
for freq in range (0,num_steps):
405
- pressure = get_pressure (node_number, freq)['x']
400
+ pressure = get_pressure_at (node_number, freq)['x']
406
401
abs_pressure = abs(pressure)
407
402
408
403
dic[start_freq] = abs_pressure
409
404
start_freq += step_val
410
405
411
- Sort the results:
406
+ Sort the results.
412
407
413
408
.. code :: python3
414
409
@@ -418,15 +413,15 @@ Sort the results:
418
413
419
414
420
415
421
- Plot the frequency response function for the selected node:
416
+ Plot the frequency response function for the selected node.
422
417
423
418
.. code :: python3
424
419
425
420
plt.plot(frf_x, frf_y, linewidth= 3.0, color='b')
426
421
427
422
# Plot the natural frequency as vertical lines on the FRF graph
428
423
for itr in range(0,6):
429
- plt.axvline(x=ev [itr], ymin=0,ymax=2, color='r', linestyle='dotted', linewidth=1)
424
+ plt.axvline(x=eigenvalues [itr], ymin=0,ymax=2, color='r', linestyle='dotted', linewidth=1)
430
425
431
426
# Name the graph and the x-axis and y-axis
432
427
plt.title("Frequency Response Function")
0 commit comments