error matrix and force vector

[mecaFEA]

🤓 Before submitting the issue

• [X] I have visited the Troubleshooting section.
• [X] I have searched among the existing issues.
• [X] I am using a Python virtual environment.
• [X] I am using the latest version of PyMAPDL and its dependencies (fully updated virtual environment). You can update them using pip install --upgrade --upgrade-strategy eager ansys-mapdl-core in your activated virtual environment.

🔍 Description of the bug

Hello,

I want to export the matrix vector and the force. Then i want to solve the linear equation by myself K*X =F.

This is the exemple:

The force vector will be like this F = [ 0 0 NODE 4 Fz ....], we will not have the node at the boundary conditions as they will embed. as there is a Force at node 4, 5 and 6 and node 1, 2, 3 are embed. The begin of my vector force normaly will be like this:

F = [ 0 0 Fz 0 0 Fz 0 0 Fz ....]

But when i import the full file with pymath library the begin of the force vector is like this :

0.000000000000000000e+00 0.000000000000000000e+00 -1.000000000000000000e+01 0.000000000000000000e+00 0.000000000000000000e+00 -1.000000000000000000e+01 0.000000000000000000e+00 0.000000000000000000e+00 0.000000000000000000e+00 0.000000000000000000e+00

The force at node are not at the good place for the others vector. So i don't know if the matrice K is well assembled

Consequently, when i solve U = K-1 * F i have not the good result. I transmit you my code that i used toi generate the mesh and solve the equation.

Thanks in advance

🕵️ Steps To Reproduce

------------------generate mesh-----------------------------

from ansys.mapdl import core as pymapdl
import numpy as np
import ansys.math.core.math as pymath
from scipy import linalg as LA
from scipy.sparse.linalg import spsolve

mapdl = pymapdl.launch_mapdl()
mapdl.prep7()
vnum0 = mapdl.block(0, 1, 0, 1, 0, 1)
mapdl.mp("EX", 1, 200e9)  # Elastic moduli in Pa (kg/(m*s**2))
mapdl.mp("DENS", 1, 7800)  # Density in kg/m3
mapdl.mp("NUXY", 1, 0.3)  # Poissons Ratio
mapdl.esize(0.5)
mapdl.et(1, 185)
mapdl.mat(1)
mapdl.type(1)
mapdl.vmesh(vnum0)
print("------------------force_node------------------")
mapdl.nsel("S","loc","z",-0.1,0.1)
mapdl.nsel("U","loc","x",-0.2,0.2)
mapdl.nplot()
mapdl.cm("force_node",'NODE')
mapdl.f(node = "all", lab = "FZ", value = -10)
mapdl.allsel()
mapdl.nsel("s","loc","x",0,0)
mapdl.cm("BC",'NODE')
BC = mapdl.mesh.nnum
mapdl.d("all", "ALL")
print("------------------BC------------------")
print(BC)
mapdl.allsel
mapdl.allsel
mapdl.eplot()
print(mapdl.slashsolu())
print(mapdl.antype(antype="MODAL"))
print(mapdl.modopt(method="LANB", nmode="10", freqb="1."))
print(mapdl.wrfull(ldstep="1"))
# # store the output of the solve command
output = mapdl.solve()

np.savetxt('BC.csv', BC, delimiter=',')

-----------------------------------solve K*U = F -------------------------

from ansys.mapdl import core as pymapdl
import numpy as np
import ansys.math.core.math as pymath
from scipy import linalg as LA
from scipy.linalg import solve

mm = pymath.AnsMath()

k = mm.stiff(fname="C:/Users/PC Prêt LAMIH/Desktop/DOCTORAT/solver_code/testt/file.full")
k = np.array(k)
k = k + k.T - np.diag(k.diagonal())
m = mm.mass(fname="C:/Users/PC Prêt LAMIH/Desktop/DOCTORAT/solver_code/testt/file.full")
m = np.array(m)
m = m + m.T - np.diag(m.diagonal())
b = mm.rhs(fname="C:/Users/PC Prêt LAMIH/Desktop/DOCTORAT/solver_code/testt/file.full", name="B")
print(b)
np.savetxt('matrix_K.csv', k, delimiter=',')
np.savetxt('matrix_M.csv', m, delimiter=',')
np.savetxt('force_vect.csv', b, delimiter=',')

U = solve(k, b)
np.savetxt('U.csv', U, delimiter=',')

💻 Which Operating System are you using?

None

🐍 Which Python version are you using?

None

📝 PyMAPDL Report

Show the Report!

```text # PASTE HERE THE OUTPUT OF `python -c "from ansys.mapdl import core as pymapdl; print(pymapdl.Report())"` here ```

📝 Installed packages

Show the installed packages!

```text # PASTE HERE THE OUTPUT OF `python -m pip freeze` here ```

📝 Logger output file

Show the logger output file.

```text # PASTE HERE THE CONTENT OF THE LOGGER OUTPUT FILE. ```

[germa89]

Assigning this to @clatapie

[mecaFEA]

I have understood that the stiffness matrix is organized arcording to the element number. In my case is element number 1 and after 6. How can i know the order in my matrix?

[mikerife]

@mecaFEA The stiffness matrix is ordered by degree-of-freedom numbering. The node and element ordering during solution can change relative to the FEM ordering. If you want to compare the MAPDL solution to the Pymath you will need to get the mapping vector and apply it to the Pymath result. See here

This is not a bug, just a slight misunderstanding. mike

[germa89]

I remember there was an example about mapping between both somewhere... @clatapie ?

[mikerife]

@clatapie @germa89 There are examples here and here (1st and 2nd examples). The 'node-to-solve' mapping matrix cannot be exported. As far as I can tell you will need to work with the Mechanical APDL dev's to be able to allow this to be exported. Or create a routine to read them from the full file. Mike

[germa89]

Closing issue then. @mecaFEA feel free to reopen if more info is required.