ansys / pymapdl

Pythonic interface to MAPDL
https://mapdl.docs.pyansys.com
MIT License
428 stars 120 forks source link

Krylov macros not reading full files in DMP #2326

Closed germa89 closed 1 year ago

germa89 commented 1 year ago

Discussed in https://github.com/ansys/pymapdl/discussions/2306

Originally posted by **dhabar95** September 5, 2023 When running the code (https://mapdl.docs.pyansys.com/version/stable/examples/extended_examples/Krylov/krylov_example.html): ````py mapdl.run('/SOLU') mapdl.antype('HARMIC') # Set options for harmonic analysis mapdl.hropt('KRYLOV') mapdl.eqslv('SPARSE') mapdl.harfrq(0,1000) # Set beginning and ending frequency mapdl.nsubst(100) # Set the number of frequency increments mapdl.wrfull(1) # Generate FULL file and stop mapdl.solve() mapdl.finish() dd = mapdl.krylov # Initialize Krylov class object Qz = dd.gensubspace(10, 500, check_orthogonality=True) print(Qz.shape) Yz = dd.solve(0, 1000, 100, ramped_load=True) print(Yz.shape) result = dd.expand(residual_computation=True, residual_algorithm="l2", return_solution = True) ```` `dd.gensubspace` is unable to locate the `file.full` ( as i am using `nproc = 4` and create `fileN.full` files). Solution tried and failed: 1. Using `COMBINE `- doesn't work for `hropt("FULL")` 2. Using `DMPOPTION `before `mapdl.solve()` - doesn't combine the `fileN.full` files 3. Using `nproc = 1` - Takes ages to solve Could you please help how to solve this problem about `nproc` and `file.full` problem?
germa89 commented 1 year ago

Running an instance (in docker) with the following command:

ansys -grpc -dmp -np 4

and then running the following code:

import os
import numpy as np
import math
import matplotlib.pyplot as plt
from ansys.mapdl.core import launch_mapdl

mapdl = launch_mapdl(start_instance=False, port=50042)
mapdl.clear()
mm = mapdl.math

# Constants
pi = np.arccos(-1)
c0 = 340                         # speed of Sound (m/s)

# Materials
rho = 1.2                        # density
c0 = 340                         # speed of Sound
frq = 1000                       # excitation freq   Hz
visco = 0.9                      # viscosity

TP = 1/frq
WL = c0 * TP
no_wl = 3                        # no of wavelengths in space

cyl_L = no_wl * WL               # length of duct
cyl_r = 0.025 * cyl_L            # cross section of duct

nelem_wl = 10                    # no of elements per wavelength
tol_elem = nelem_wl * no_wl      # total number of elements across length

mapdl.prep7()

mapdl.et(1,'FLUID220', kop2=1)   # uncoupled acoustic element without FSIs
mapdl.mp("DENS", 1, rho)
mapdl.mp("SONC", 1, c0)
mapdl.mp("VISC", 1, visco)

# Set back to default
mapdl.csys(0)

# Rotate working plane for the cylinder generation
mapdl.wpcsys(-1)
mapdl.wprota(thzx=90)

# Generate a circular area with a specified radius
mapdl.cyl4(0, 0, cyl_r)

mapdl.wpcsys(-1)

# Extrude the circular area to generate a cylinder of specified length
mapdl.vext("ALL", dx=cyl_L)

# Split the cylinder into four segments to create a more uniform mesh
mapdl.vsbw("ALL", keep='DELETE')
mapdl.wprota(thzx=90)
mapdl.vsbw("ALL", keep='DELETE')

mapdl.wpcsys(-1)

# Create a component with the created volume
mapdl.cm('cm1', 'volu')

# Select material and type
mapdl.mat(1)
mapdl.type(1)

# Select volume to mesh
mapdl.cmsel("S", "cm1")

# Select lines belonging to the volume
mapdl.aslv()
mapdl.lsla()

# Unselect lines at the top and bottom faces
mapdl.lsel("U", 'loc', 'x', 0)
mapdl.lsel("U", 'loc', 'x', cyl_L)

# Apply length constraint
mapdl.lesize('ALL',ndiv = tol_elem)
mapdl.lsla()

# Mesh
mapdl.vsweep('ALL')
mapdl.allsel()

# Plot the FE model
mapdl.eplot()

# Select areas to apply pressure to
mapdl.cmsel("S", "cm1")
mapdl.aslv()
mapdl.asel('R',"EXT")  # select external areas
mapdl.asel('R',"LOC","x",0)
mapdl.nsla('S',1)

# Apply pressure
mapdl.d('ALL','PRES', 1)

# Select nodes on the areas where impedance is to be applied
mapdl.cmsel("S", "cm1")
mapdl.aslv()
mapdl.asel('R',"EXT")
mapdl.asel('R',"LOC","x",cyl_L)
mapdl.nsla("S",1)

# Apply impedance
mapdl.sf("ALL","IMPD",1000)
mapdl.allsel()

# Modal Analysis
mapdl.slashsolu()
nev = 10 # Get the first 10 modes
output = mapdl.modal_analysis("DAMP", nmode=nev)
mapdl.finish()
mm.free()

k = mm.stiff(fname=f"{mapdl.jobname}.full")
M = mm.mass(fname=f"{mapdl.jobname}.full")
A = mm.mat(k.nrow, nev)
eigenvalues = mm.eigs(nev, k, M, phi=A, fmin=1.0)

for each_freq in range(10):
     print(f"Freq = {eigenvalues[each_freq]:8.2f} Hz") # Eigenfrequency (Hz)

mapdl.antype('HARMIC')  # Set options for harmonic analysis
mapdl.hropt('KRYLOV')
mapdl.eqslv('SPARSE')
mapdl.harfrq(0,1000)    # Set beginning and ending frequency
mapdl.nsubst(100)       # Set the number of frequency increments
mapdl.wrfull(1)         # Generate FULL file and stop
mapdl.solve()
mapdl.finish()

dd = mapdl.krylov       # Initialize Krylov class object

Qz = dd.gensubspace(10, 500, check_orthogonality=True)

print(Qz.shape)

Did NOT raise any issue.

dhabar95 commented 1 year ago

Thanks a lot for the solution! Helped a lot :smiley:

germa89 commented 1 year ago

Let's close this issue then.