ansys / pymapdl

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

Matrix multiplication from MAPDL doesn't match matrix multiplication from scipy and numpy #1157

Closed raphael-31 closed 10 months ago

raphael-31 commented 2 years ago

Describe the bug Hi, I am working on the resolution of eigen problem to find the added mass of a system. I have a problem with a simple matrix multiplication. The use of the .dot method from the math module of MAPDL dosen't give the same result of multplying matrix with numpy. So I canno't verify that the eigenvalue is verified with the APDL math matrix converted as array.

To Reproduce Steps to reproduce the behavior:

  1. Load the stiff and the mass matrix from an ansys APDL model
  2. Solve the eigen problem with the eig method to obtain the eigen vectors phis
  3. Choose one vector phi
  4. Multiply : K.dot(phi)
  5. Convert : K and phi into array, multiply K and phi with numpy or scipy methods
  6. The multiplication give not the same vector, so the eigenproblem is not verified with numpy matrix

Expected behavior It should be expected that the multiplication gives the same result and the eigen problem is verified whether the matrix are under the mapdl form or the numpy form

Screenshots The multiplication with ADPL: image

The multiplication with numpy: image

Only the first component of the vectors are equal the others component are not in the same order of magnitude

System Information:

My code:

Toggle me! ```Python print("version of ansys.mapdl") mapdl = launch_mapdl(run_location = folder) print(mapdl) print("----" + "case : " + args.test_case + "----") folder = root / folder_name file_name_full = "cylinder_no_water.full" data_file_name = folder / "cylinder_no_water.dat" folder = str(folder) mm = mapdl.math mapdl.finish() mm.free() K = mm.stiff(fname=file_name_full) M = mm.mass(fname=file_name_full) n_ev = 5 eigen_vectors = mm.zeros(K.nrow, n_ev) ev = mm.eigs(n_ev, K, M, phi=eigen_vectors, fmin = 1) index_ev = 2 eigen_vec = eigen_vectors[index_ev] freq = ev[index_ev] omega = 2 * np.pi * freq lam = omega * omega #-----------APDL MATRIX-------------- print("using apdl matrix") Kphi = K.dot(eigen_vec) Mphi = M.dot(eigen_vec) Mphi *= lam res = (Kphi - Mphi).norm() res = res / Kphi.norm() print(Kphi) print(res) #-----------scipy MATRIX-------------- print("using scipy matrix") eigen_vec_sci = eigen_vec.asarray() K_sci = K.asarray() M_sci = M.asarray() Kphi = K_sci * eigen_vec_sci.T print(Kphi) #-----------numpy MATRIX-------------- print("using numpy matrix") eigen_vec_np = eigen_vec.asarray() K_np = K.todense() M_np = M.todense() Kphi = np.matmul(K_np, eigen_vec_np) Mphi = np.matmul(M_np, eigen_vec_np) Mphi *= lam res = np.linalg.norm(Kphi - Mphi) res = res / np.linalg.norm(Kphi) print(res) mapdl.exit() ``` ``` version of ansys.mapdl Product: Ansys Mechanical Enterprise Academic Student MAPDL Version: 22.1 ansys.mapdl Version: 0.60.6 ----case : test-extract-matrix-basic---- using apdl matrix OSBOEH : Size : 880 -8.534e+04 5.069e+05 -1.244e+05 4.607e+05 -2.335e+05 < 5 9.160e+05 -1.717e+05 1.003e+06 -2.112e+05 9.075e+05 < 10 -1.730e+05 9.864e+05 -1.693e+05 8.879e+05 -8.546e+04 < 15 6.645e+05 -5.280e+04 6.680e+05 -1.495e+05 8.769e+05 < 20 -7.429e+04 8.255e+05 6.897e+04 6.731e+05 4.436e+04 < 25 3.364e+05 -3.095e+04 4.079e+05 -9.283e+04 8.395e+05 < 30 3.929e+04 6.730e+05 -1.114e+05 8.527e+05 9.292e+03 < 35 6.722e+05 -1.303e+05 8.652e+05 -2.131e+04 6.706e+05 < 40 -1.898e+05 8.981e+05 -1.196e+05 6.597e+05 -1.553e+05 < 45 6.538e+05 -1.931e+05 6.464e+05 -1.096e+05 3.204e+05 < 50 9.978e+04 3.535e+05 1.340e+05 3.733e+05 6.374e+04 < 55 3.324e+05 -1.488e+04 2.848e+05 2.562e+04 3.096e+05 < 60 -2.159e+04 1.337e+05 2.595e+05 4.418e+05 1.398e+05 < 65 2.261e+05 2.292e+05 4.258e+05 1.985e+05 4.092e+05 < 70 1.668e+05 3.917e+05 4.119e+05 1.048e+05 4.341e+05 < 75 1.333e+05 3.637e+05 4.455e+04 3.885e+05 7.531e+04 < 80 2.778e+05 -5.792e+04 3.085e+05 -2.180e+04 1.280e+05 < 85 -4.157e+04 3.371e+05 1.227e+04 4.760e+05 1.884e+05 < 90 2.447e+05 1.033e+05 4.554e+05 1.611e+05 6.750e+05 < 95 -3.600e+04 3.389e+05 -8.221e+03 6.702e+05 -6.532e+04 < 100 6.646e+05 -9.484e+04 6.582e+05 -1.248e+05 6.507e+05 < 105 -1.556e+05 6.420e+05 -1.873e+05 6.320e+05 -2.202e+05 < 110 6.204e+05 -2.547e+05 6.072e+05 -2.908e+05 2.987e+05 < 115 -1.579e+05 4.352e+05 -1.947e+05 8.677e+05 -3.733e+05 < 120 8.628e+05 -3.499e+05 8.424e+05 -2.841e+05 8.339e+05 < 125 -2.632e+05 8.501e+05 -3.054e+05 8.569e+05 -3.273e+05 < 130 8.036e+05 -2.015e+05 3.980e+05 -9.386e+04 8.145e+05 < 135 -2.221e+05 8.246e+05 -2.426e+05 4.874e+05 -1.635e+05 < 140 9.640e+05 -3.263e+05 8.765e+05 -3.154e+05 8.575e+05 < 145 -3.123e+05 9.474e+05 -3.248e+05 8.378e+05 -3.088e+05 < 150 9.129e+05 -3.208e+05 8.950e+05 -3.183e+05 9.304e+05 < 155 -3.230e+05 8.175e+05 -3.050e+05 4.018e+05 -1.512e+05 < 160 6.973e+05 -3.656e+05 3.418e+05 -1.878e+05 7.175e+05 < 165 -3.504e+05 7.371e+05 -3.347e+05 7.562e+05 -3.182e+05 < 170 7.750e+05 -3.009e+05 7.934e+05 -2.824e+05 8.115e+05 < 175 -2.627e+05 8.293e+05 -2.415e+05 8.470e+05 -2.186e+05 < 180 4.293e+05 -1.012e+05 4.984e+05 -2.874e+05 4.887e+05 < 185 -3.201e+05 5.256e+05 -1.777e+05 5.168e+05 -2.165e+05 < 190 5.341e+05 -1.362e+05 2.698e+05 -5.356e+04 4.787e+05 < 195 -3.516e+05 5.077e+05 -2.529e+05 4.575e+05 -4.117e+05 < 200 2.251e+05 -2.156e+05 4.683e+05 -3.819e+05 1.439e+05 < 205 -4.648e+05 7.318e+04 -2.438e+05 1.401e+05 -4.302e+05 < 210 1.318e+05 -3.583e+05 1.360e+05 -3.948e+05 1.272e+05 < 215 -3.202e+05 1.223e+05 -2.802e+05 1.170e+05 -2.378e+05 < 220 1.049e+05 -1.447e+05 1.113e+05 -1.928e+05 5.022e+04 < 225 -5.548e+04 -1.796e+05 -5.278e+05 -8.392e+04 -2.732e+05 < 230 -1.970e+05 -4.996e+05 -2.144e+05 -4.706e+05 -2.318e+05 < 235 -4.404e+05 -2.495e+05 -4.088e+05 -2.676e+05 -3.754e+05 < 240 -2.861e+05 -3.398e+05 -3.052e+05 -3.018e+05 -3.250e+05 < 245 -2.609e+05 -1.693e+05 -1.161e+05 -4.736e+05 -5.683e+05 < 250 -4.493e+05 -5.814e+05 -4.975e+05 -5.545e+05 -2.165e+05 < 255 -2.949e+05 -5.214e+05 -5.398e+05 -5.452e+05 -5.239e+05 < 260 -5.693e+05 -5.066e+05 -6.185e+05 -4.672e+05 -5.937e+05 < 265 -4.878e+05 -6.438e+05 -4.447e+05 -3.305e+05 -2.143e+05 < 270 -6.421e+05 -5.964e+05 -6.211e+05 -5.912e+05 -6.626e+05 < 275 -6.009e+05 -3.034e+05 -2.937e+05 -6.827e+05 -6.049e+05 < 280 -7.026e+05 -6.081e+05 -7.222e+05 -6.107e+05 -7.418e+05 < 285 -6.125e+05 -7.613e+05 -6.135e+05 -7.808e+05 -6.135e+05 < 290 -3.969e+05 -3.066e+05 -6.807e+05 -5.219e+05 -3.374e+05 < 295 -2.542e+05 -6.892e+05 -5.418e+05 -6.969e+05 -5.613e+05 < 300 -7.040e+05 -5.805e+05 -7.104e+05 -5.995e+05 -7.161e+05 < 305 -6.185e+05 -7.210e+05 -6.375e+05 -7.251e+05 -6.565e+05 < 310 -7.282e+05 -6.757e+05 -3.649e+05 -3.443e+05 -6.446e+05 < 315 -3.541e+05 -3.252e+05 -1.683e+05 -6.356e+05 -3.801e+05 < 320 -6.257e+05 -4.059e+05 -6.149e+05 -4.317e+05 -6.030e+05 < 325 -4.578e+05 -5.898e+05 -4.842e+05 -5.405e+05 -5.675e+05 < 330 -5.587e+05 -5.389e+05 -5.751e+05 -5.112e+05 -2.636e+05 < 335 -2.935e+05 -5.497e+05 -9.625e+04 -2.831e+05 -4.088e+04 < 340 -5.246e+05 -1.179e+05 -4.987e+05 -1.396e+05 -4.717e+05 < 345 -1.616e+05 -4.433e+05 -1.841e+05 -4.132e+05 -2.071e+05 < 350 -3.810e+05 -2.310e+05 -3.464e+05 -2.558e+05 -3.092e+05 < 355 -2.818e+05 -1.415e+05 -1.499e+05 -4.371e+05 2.133e+05 < 360 -2.296e+05 1.097e+05 -4.035e+05 2.041e+05 -3.693e+05 < 365 1.946e+05 -3.338e+05 1.846e+05 -2.969e+05 1.741e+05 < 370 -2.582e+05 1.630e+05 -2.172e+05 1.511e+05 -1.736e+05 < 375 1.383e+05 -1.271e+05 1.245e+05 -4.724e+04 5.739e+04 < 380 -3.357e+05 5.150e+05 -1.781e+05 2.554e+05 -3.046e+05 < 385 5.210e+05 -2.730e+05 5.264e+05 -2.403e+05 5.314e+05 < 390 -2.064e+05 5.358e+05 -1.709e+05 5.396e+05 -1.334e+05 < 395 5.428e+05 -9.369e+04 5.454e+05 -5.138e+04 5.473e+05 < 400 -1.088e+04 2.741e+05 -2.528e+05 7.451e+05 -1.324e+05 < 405 3.665e+05 -2.345e+05 7.626e+05 -2.159e+05 7.795e+05 < 410 -1.966e+05 7.958e+05 -1.765e+05 8.116e+05 -1.553e+05 < 415 8.268e+05 -1.329e+05 8.416e+05 -1.091e+05 8.559e+05 < 420 -8.364e+04 8.698e+05 -3.289e+04 4.394e+05 -1.742e+05 < 425 8.549e+05 -8.691e+04 4.203e+05 -1.747e+05 8.756e+05 < 430 -1.750e+05 8.955e+05 -1.751e+05 9.148e+05 -1.750e+05 < 435 9.335e+05 -1.746e+05 9.516e+05 -1.739e+05 9.692e+05 < 440 -1.730e+05 9.864e+05 -1.717e+05 1.003e+06 -8.534e+04 < 445 5.069e+05 -1.114e+05 8.527e+05 -9.283e+04 8.395e+05 < 450 -1.303e+05 8.652e+05 -1.495e+05 8.769e+05 -1.693e+05 < 455 8.879e+05 -1.898e+05 8.981e+05 -2.112e+05 9.075e+05 < 460 -2.335e+05 9.160e+05 -1.244e+05 4.607e+05 4.575e+05 < 465 -4.117e+05 1.439e+05 -4.648e+05 7.318e+04 -2.438e+05 < 470 2.251e+05 -2.156e+05 4.683e+05 -3.819e+05 1.401e+05 < 475 -4.302e+05 4.787e+05 -3.516e+05 1.360e+05 -3.948e+05 < 480 4.887e+05 -3.201e+05 1.318e+05 -3.583e+05 4.984e+05 < 485 -2.874e+05 1.272e+05 -3.202e+05 5.077e+05 -2.529e+05 < 490 1.223e+05 -2.802e+05 5.168e+05 -2.165e+05 1.170e+05 < 495 -2.378e+05 5.256e+05 -1.777e+05 1.113e+05 -1.928e+05 < 500 5.341e+05 -1.362e+05 1.049e+05 -1.447e+05 2.698e+05 < 505 -5.356e+04 5.022e+04 -5.548e+04 -2.495e+05 -4.088e+05 < 510 -2.318e+05 -4.404e+05 -3.052e+05 -3.018e+05 -2.861e+05 < 515 -3.398e+05 -3.250e+05 -2.609e+05 -1.796e+05 -5.278e+05 < 520 -8.392e+04 -2.732e+05 -1.970e+05 -4.996e+05 -2.144e+05 < 525 -4.706e+05 -2.676e+05 -3.754e+05 -1.693e+05 -1.161e+05 < 530 -4.493e+05 -5.814e+05 -2.165e+05 -2.949e+05 -4.736e+05 < 535 -5.683e+05 -4.975e+05 -5.545e+05 -5.214e+05 -5.398e+05 < 540 -5.452e+05 -5.239e+05 -5.693e+05 -5.066e+05 -5.937e+05 < 545 -4.878e+05 -6.185e+05 -4.672e+05 -6.438e+05 -4.447e+05 < 550 -3.305e+05 -2.143e+05 -6.211e+05 -5.912e+05 -3.034e+05 < 555 -2.937e+05 -6.421e+05 -5.964e+05 -6.626e+05 -6.009e+05 < 560 -6.827e+05 -6.049e+05 -7.026e+05 -6.081e+05 -7.418e+05 < 565 -6.125e+05 -7.222e+05 -6.107e+05 -7.613e+05 -6.135e+05 < 570 -7.808e+05 -6.135e+05 -3.969e+05 -3.066e+05 -6.807e+05 < 575 -5.219e+05 -3.374e+05 -2.542e+05 -6.892e+05 -5.418e+05 < 580 -6.969e+05 -5.613e+05 -7.040e+05 -5.805e+05 -7.104e+05 < 585 -5.995e+05 -7.161e+05 -6.185e+05 -7.210e+05 -6.375e+05 < 590 -7.251e+05 -6.565e+05 -7.282e+05 -6.757e+05 -3.649e+05 < 595 -3.443e+05 -6.446e+05 -3.541e+05 -3.252e+05 -1.683e+05 < 600 -6.356e+05 -3.801e+05 -6.149e+05 -4.317e+05 -6.257e+05 < 605 -4.059e+05 -6.030e+05 -4.578e+05 -5.898e+05 -4.842e+05 < 610 -5.751e+05 -5.112e+05 -5.587e+05 -5.389e+05 -5.405e+05 < 615 -5.675e+05 -2.636e+05 -2.935e+05 -5.497e+05 -9.625e+04 < 620 -2.831e+05 -4.088e+04 -5.246e+05 -1.179e+05 -4.987e+05 < 625 -1.396e+05 -4.717e+05 -1.616e+05 -4.433e+05 -1.841e+05 < 630 -4.132e+05 -2.071e+05 -3.810e+05 -2.310e+05 -3.464e+05 < 635 -2.558e+05 -3.092e+05 -2.818e+05 -1.415e+05 -1.499e+05 < 640 -4.371e+05 2.133e+05 -2.296e+05 1.097e+05 -4.035e+05 < 645 2.041e+05 -3.693e+05 1.946e+05 -3.338e+05 1.846e+05 < 650 -2.969e+05 1.741e+05 -2.582e+05 1.630e+05 -2.172e+05 < 655 1.511e+05 -1.736e+05 1.383e+05 -1.271e+05 1.245e+05 < 660 -4.724e+04 5.739e+04 -3.357e+05 5.150e+05 -1.781e+05 < 665 2.554e+05 -3.046e+05 5.210e+05 -2.730e+05 5.264e+05 < 670 -2.403e+05 5.314e+05 -2.064e+05 5.358e+05 -1.709e+05 < 675 5.396e+05 -1.334e+05 5.428e+05 -9.369e+04 5.454e+05 < 680 -5.138e+04 5.473e+05 -1.088e+04 2.741e+05 -2.528e+05 < 685 7.451e+05 -1.324e+05 3.665e+05 -2.345e+05 7.626e+05 < 690 -2.159e+05 7.795e+05 -1.966e+05 7.958e+05 -1.765e+05 < 695 8.116e+05 -1.553e+05 8.268e+05 -1.329e+05 8.416e+05 < 700 -1.091e+05 8.559e+05 -8.364e+04 8.698e+05 -3.289e+04 < 705 4.394e+05 -1.746e+05 9.516e+05 -1.750e+05 9.335e+05 < 710 -1.739e+05 9.692e+05 -1.742e+05 8.549e+05 -8.691e+04 < 715 4.203e+05 -1.747e+05 8.756e+05 -1.750e+05 8.955e+05 < 720 -1.751e+05 9.148e+05 -7.429e+04 8.255e+05 -3.095e+04 < 725 4.079e+05 6.897e+04 6.731e+05 4.436e+04 3.364e+05 < 730 3.929e+04 6.730e+05 9.292e+03 6.722e+05 -2.131e+04 < 735 6.706e+05 -8.546e+04 6.645e+05 -5.280e+04 6.680e+05 < 740 -1.196e+05 6.597e+05 -1.553e+05 6.538e+05 -1.931e+05 < 745 6.464e+05 -1.096e+05 3.204e+05 2.292e+05 4.258e+05 < 750 2.595e+05 4.418e+05 1.985e+05 4.092e+05 1.398e+05 < 755 2.261e+05 1.668e+05 3.917e+05 1.340e+05 3.733e+05 < 760 9.978e+04 3.535e+05 6.374e+04 3.324e+05 2.562e+04 < 765 3.096e+05 -1.488e+04 2.848e+05 -2.159e+04 1.337e+05 < 770 4.760e+05 1.884e+05 2.447e+05 1.033e+05 4.554e+05 < 775 1.611e+05 4.341e+05 1.333e+05 4.119e+05 1.048e+05 < 780 3.885e+05 7.531e+04 3.637e+05 4.455e+04 3.371e+05 < 785 1.227e+04 3.085e+05 -2.180e+04 2.778e+05 -5.792e+04 < 790 1.280e+05 -4.157e+04 6.750e+05 -3.600e+04 3.389e+05 < 795 -8.221e+03 6.702e+05 -6.532e+04 6.646e+05 -9.484e+04 < 800 6.582e+05 -1.248e+05 6.507e+05 -1.556e+05 6.420e+05 < 805 -1.873e+05 6.320e+05 -2.202e+05 6.204e+05 -2.547e+05 < 810 6.072e+05 -2.908e+05 2.987e+05 -1.579e+05 8.036e+05 < 815 -2.015e+05 3.980e+05 -9.386e+04 8.145e+05 -2.221e+05 < 820 8.246e+05 -2.426e+05 8.339e+05 -2.632e+05 8.424e+05 < 825 -2.841e+05 8.501e+05 -3.054e+05 8.569e+05 -3.273e+05 < 830 8.628e+05 -3.499e+05 8.677e+05 -3.733e+05 4.352e+05 < 835 -1.947e+05 8.175e+05 -3.050e+05 4.018e+05 -1.512e+05 < 840 8.378e+05 -3.088e+05 8.575e+05 -3.123e+05 8.765e+05 < 845 -3.154e+05 8.950e+05 -3.183e+05 9.129e+05 -3.208e+05 < 850 9.304e+05 -3.230e+05 9.474e+05 -3.248e+05 9.640e+05 < 855 -3.263e+05 4.874e+05 -1.635e+05 6.973e+05 -3.656e+05 < 860 3.418e+05 -1.878e+05 7.175e+05 -3.504e+05 7.371e+05 < 865 -3.347e+05 7.562e+05 -3.182e+05 7.750e+05 -3.009e+05 < 870 7.934e+05 -2.824e+05 8.115e+05 -2.627e+05 8.293e+05 < 875 -2.415e+05 8.470e+05 -2.186e+05 4.293e+05 -1.012e+05 < 880 2.4333951633804383e-08 using scipy matrix [-8.53352496e+04 -1.74778086e+08 4.29756066e+08 -3.07402985e+09 -2.99527961e+09 1.04422745e+10 -1.27468682e+09 7.12713273e+09 -2.89098266e+09 1.03605155e+10 -1.29239990e+09 6.98674370e+09 -1.69336506e+05 -6.90814908e+08 1.24915996e+09 -4.79871822e+09 -2.60856496e+09 7.55466746e+09 -2.18134080e+09 5.95297944e+09 -7.42874076e+04 -3.07963450e+08 9.01959712e+08 -3.55050182e+09 -1.79433405e+09 7.58967708e+09 -2.09964573e+09 7.18606481e+09 -1.95077262e+09 8.92997853e+09 -1.26103040e+09 3.38942185e+09 -2.09396524e+09 9.02864885e+09 -1.36232876e+09 3.27451375e+09 -4.99931445e+09 1.86425037e+10 -3.90377565e+09 1.09660494e+10 -5.45726806e+09 2.03714496e+10 -2.24077780e+09 5.07972595e+09 -2.38207400e+09 4.94026068e+09 -2.52736520e+09 4.78442452e+09 -2.92340468e+09 5.10218615e+09 3.25210432e+07 5.82904198e+08 -6.23997351e+08 3.00198890e+09 -6.52878377e+08 2.58183605e+09 9.61832729e+06 -1.73515517e+08 -1.82608980e+09 4.36367508e+09 -1.17615408e+09 1.76024302e+09 7.95292453e+07 1.63898583e+09 -2.72876723e+06 3.40289695e+09 -5.10246838e+07 3.79633963e+09 -1.68535125e+08 3.57535400e+09 -8.35816269e+08 5.56886401e+09 7.43195355e+07 2.69070671e+09 2.22680507e+09 2.07359212e+09 5.46861699e+07 2.35301571e+09 4.30161381e+09 7.50550213e+08 2.82782511e+07 1.75545662e+09 1.90608733e+09 6.27287077e+08 1.93748193e+09 -6.21796477e+08 4.07377273e+09 -5.78034603e+07 1.09094738e+08 3.14168906e+09 2.22323523e+09 1.25550783e+09 4.57684635e+09 1.85058673e+09 1.22185104e+08 4.23982378e+09 4.79170233e+09 -2.90645623e+08 5.13349933e+09 1.80201523e+09 5.18855062e+09 1.62017225e+09 5.23928852e+09 1.43415355e+09 5.28554663e+09 1.24194318e+09 5.32708399e+09 1.04158943e+09 5.36359128e+09 8.31188043e+08 5.39469593e+09 6.08867942e+08 5.41996648e+09 3.72779010e+08 5.53682952e+09 -1.58775354e+09 1.81887743e+08 2.59428347e+09 8.28121949e+09 5.50873143e+08 8.17547925e+09 6.44343191e+08 9.15263236e+07 4.49673520e+09 7.70140410e+09 9.58844485e+08 7.84441831e+09 9.74506271e+08 1.59549646e+10 -2.94979015e+09 1.19280159e+08 4.29218626e+09 6.87927676e+09 -1.13764953e+09 7.28697148e+09 1.26359916e+09 1.49011061e+10 -2.18126106e+09 1.53982561e+08 2.25611442e+09 9.48817308e+09 6.72143993e+08 9.33245977e+07 3.39055728e+09 8.47640785e+09 5.09875533e+08 9.32970060e+09 6.43836722e+08 8.29194298e+09 4.83531221e+08 8.69347502e+07 3.53859432e+09 1.73846916e+10 -2.29982182e+09 1.80823896e+10 -2.34477939e+09 8.10172272e+09 4.56924009e+08 7.71416119e+09 -1.44236268e+09 8.13258728e+07 1.36796352e+09 6.86845267e+09 -1.45572001e+09 7.08348252e+09 -3.94753202e+08 7.26742884e+09 -3.25792247e+08 7.44527299e+09 -2.55591868e+08 7.61787726e+09 -1.83472976e+08 7.78593188e+09 -1.08776901e+08 7.94997682e+09 -3.08602620e+07 8.11042006e+09 5.09094542e+07 8.26755300e+09 1.37157725e+08 8.43166906e+09 -3.07175977e+08 5.53945599e+07 -5.94534698e+07 4.79296371e+09 -8.22767484e+08 5.29991569e+07 -3.80308379e+07 5.07365578e+09 -4.89768493e+08 5.16815596e+09 -4.84878152e+08 5.20970948e+09 -9.35040346e+07 4.69385315e+09 -9.23080715e+08 9.81639272e+09 -1.37892057e+09 5.60677497e+07 -9.26093691e+07 4.43225716e+09 -1.41015703e+09 9.00767971e+09 -2.14946264e+09 3.11174972e+07 -4.10699496e+08 1.00207259e+09 -1.48393928e+09 8.98778456e+08 -1.58960918e+09 2.84019447e+07 -3.47180128e+08 1.74865750e+09 -2.33715637e+09 8.94889023e+08 -1.21512267e+09 8.90958152e+08 -1.08219890e+09 8.85530144e+08 -9.43006201e+08 2.01668353e+07 -2.37884697e+08 1.70368733e+09 -1.05258791e+09 8.02421440e+08 -2.64552231e+08 7.15622976e+06 5.73828428e+08 -2.51435436e+09 -1.78786625e+09 -2.68143527e+09 -1.50096224e+09 -2.78871259e+09 -1.35263041e+09 -2.89414674e+09 -1.20156942e+09 -2.99838989e+09 -1.04617941e+09 -3.10201604e+09 -8.84914414e+08 -3.20553190e+09 -7.16268654e+08 -3.30938617e+09 -5.38765029e+08 -3.41397733e+09 -3.50945189e+08 -3.52492279e+09 -8.90559829e+08 -2.05387426e+07 2.43344244e+09 -5.39389368e+09 -1.33107459e+09 -5.51785329e+09 -1.11200488e+09 -5.20536034e+09 -2.31590817e+09 -5.68875963e+09 -9.75506896e+08 -5.85636184e+09 -8.34808245e+08 -6.02148330e+09 -6.88575818e+08 -4.25944451e+07 3.32285711e+09 -1.25251099e+10 -4.14689281e+09 -6.50841422e+09 -2.03730621e+08 -6.57541433e+09 -2.02247125e+09 -4.52153593e+07 3.93995279e+09 -6.53452567e+09 -1.27448490e+09 -6.68195726e+09 -1.08809364e+09 -6.46135227e+09 -2.87390923e+09 -6.85825169e+09 -1.04041168e+09 -7.02988301e+09 -9.90275015e+08 -7.19758055e+09 -9.37196839e+08 -7.36193081e+09 -8.80684821e+08 -7.52339508e+09 -8.20239646e+08 -7.68232439e+09 -7.55353696e+08 -7.63305079e+09 -3.46568777e+09 -6.82879786e+07 4.38835570e+09 -6.16670704e+09 -3.07259075e+09 -6.09074265e+09 -1.13107197e+09 -6.21894402e+09 -1.22978747e+09 -6.34159342e+09 -1.32841979e+09 -6.45913190e+09 -1.42786275e+09 -6.57186579e+09 -1.52894456e+09 -6.67998267e+09 -1.63243928e+09 -6.78356469e+09 -1.73907642e+09 -6.88259973e+09 -1.84954912e+09 -6.69225581e+09 -4.60277188e+09 -9.37776409e+07 3.95191645e+09 -4.70829271e+09 -2.46177165e+09 -4.56039225e+09 -7.93992470e+08 -4.59767517e+09 -1.02158744e+09 -4.63043302e+09 -1.25146606e+09 -4.65843749e+09 -1.48576388e+09 -4.68138797e+09 -1.72650565e+09 -6.31448128e+07 3.30999683e+09 -4.63540372e+09 -2.44176649e+09 -9.26379025e+09 -7.68255187e+09 -4.42101039e+09 -4.54571675e+09 -1.14122677e+08 2.86889713e+09 -2.80084535e+09 -7.82099856e+08 -2.66504274e+09 4.00901071e+08 -2.61028150e+09 1.43488011e+08 -2.55249512e+09 -1.18210337e+08 -2.49081706e+09 -3.86734024e+08 -2.42437763e+09 -6.64519762e+08 -2.35229997e+09 -9.53923972e+08 -2.27369646e+09 -1.25724240e+09 -2.18766569e+09 -1.57672685e+09 -1.87549609e+09 -2.64203849e+09 -1.24029156e+08 1.67849325e+09 -1.19853017e+09 1.80523257e+09 -1.16531683e+09 2.57100807e+09 -1.05267517e+09 2.41485897e+09 -9.38374074e+08 2.25331883e+09 -8.21181766e+08 2.08454921e+09 -6.99912313e+08 1.90674438e+09 -5.73414618e+08 1.71811882e+09 -4.40563062e+08 1.51689660e+09 -3.00249526e+08 1.30130212e+09 -7.82038540e+07 9.75944933e+08 -1.18945133e+08 7.25503053e+08 -4.04371633e+08 4.67927759e+09 -5.32351458e+08 5.18169474e+09 -4.16936466e+08 5.20607423e+09 -3.00663837e+08 5.22470483e+09 -1.82352761e+08 5.23713199e+09 -6.08799665e+07 5.24282753e+09 6.48323611e+07 5.24119466e+09 1.95826797e+08 5.23157204e+09 3.33119359e+08 5.21323717e+09 3.77535807e+08 5.23078206e+09 -9.68745092e+07 8.47733462e+07 -4.88785360e+08 6.97978523e+09 -7.77795344e+08 7.27355110e+09 -7.13425779e+08 7.44961878e+09 -6.48339398e+08 7.61967994e+09 -5.81812143e+08 7.78445874e+09 -5.13150944e+08 7.94451624e+09 -4.41686565e+08 8.10027029e+09 -3.66767565e+08 8.25201216e+09 -2.87755187e+08 8.39992052e+09 -4.54935609e+08 8.53623369e+09 -5.95249606e+07 -2.58831993e+08 -1.08548222e+09 7.98890184e+09 -1.47188554e+09 8.03106380e+09 -1.48958205e+09 8.23641280e+09 -1.50629376e+09 8.43503896e+09 -1.52199816e+09 8.62786996e+09 -1.53665417e+09 8.81564352e+09 -1.55020364e+09 8.99893113e+09 -1.56257254e+09 9.17815806e+09 -1.57367204e+09 9.35362012e+09 -1.91824057e+09 9.61618385e+09 -1.09193251e+07 -3.65025377e+08 -2.08800931e+09 7.45923231e+09 -2.13032005e+09 7.47042165e+09 -2.22412610e+09 7.55091382e+09 -2.31802463e+09 7.62345631e+09 -2.41256609e+09 7.68779635e+09 -2.50825329e+09 7.74356472e+09 -2.60554853e+09 7.79028619e+09 -3.03528436e+09 8.19370208e+09 4.57527425e+05 -4.11660656e+05 -2.18475346e+09 7.00934646e+08 2.10123775e+09 -1.98645882e+09 4.21983835e+09 -1.07381554e+09 4.88337230e+09 -1.81282275e+09 -3.90351745e+08 -1.01116017e+09 4.97211584e+09 -1.67385091e+09 -4.08061476e+08 -9.35519184e+08 5.05697628e+09 -1.53175172e+09 -4.26367999e+08 -8.59618464e+08 5.13816384e+09 -1.38494138e+09 -4.45358551e+08 -7.82528204e+08 5.21580069e+09 -1.23187725e+09 -4.65123932e+08 -7.03360938e+08 5.28993064e+09 -1.07104586e+09 -4.85758170e+08 -6.21261988e+08 5.36052717e+09 -9.00952764e+08 -5.07358315e+08 -5.35401406e+08 5.42750026e+09 -7.20113672e+08 -5.30024330e+08 -4.44967130e+08 5.49042966e+09 -5.27008663e+08 8.02421440e+08 -2.64552231e+08 -3.27662086e+06 9.15156584e+08 -3.03473496e+09 -1.01755149e+09 -1.22246103e+07 1.17839084e+09 -3.34370727e+09 -5.01963464e+08 -3.41397733e+09 -3.50945189e+08 7.15622973e+06 5.73828428e+08 -2.51435436e+09 -1.78786625e+09 -2.68143527e+09 -1.50096224e+09 -5.72299894e+09 -3.27306266e+09 -6.33597197e+09 -2.56858946e+09 -3.52492279e+09 -8.90559829e+08 -1.67729730e+07 2.28304321e+09 -5.20536034e+09 -2.31590817e+09 -5.34267488e+09 -1.24567535e+09 -5.51785329e+09 -1.11200488e+09 -5.68875963e+09 -9.75506896e+08 -5.85636184e+09 -8.34808245e+08 -6.02148330e+09 -6.88575818e+08 -6.18482220e+09 -5.35505811e+08 -6.34696780e+09 -3.74314552e+08 -6.50841422e+09 -2.03730621e+08 -6.57541433e+09 -2.02247125e+09 -4.19099857e+07 3.78975634e+09 -6.46135227e+09 -2.87390923e+09 -6.50010638e+09 -1.13380379e+09 -6.68195726e+09 -1.08809364e+09 -6.85825169e+09 -1.04041168e+09 -7.02988301e+09 -9.90275016e+08 -6.01330192e+07 4.65282245e+09 -1.45146388e+10 -6.44199321e+09 -7.52339508e+09 -8.20239647e+08 -7.68232439e+09 -7.55353696e+08 -7.63305079e+09 -3.46568777e+09 -6.82879785e+07 4.38835570e+09 -6.16670704e+09 -3.07259075e+09 -6.09074265e+09 -1.13107197e+09 -6.21894402e+09 -1.22978747e+09 -6.34159342e+09 -1.32841979e+09 -6.45913190e+09 -1.42786275e+09 -6.57186579e+09 -1.52894456e+09 -6.67998267e+09 -1.63243928e+09 -6.78356469e+09 -1.73907642e+09 -6.88259973e+09 -1.84954912e+09 -6.69225581e+09 -4.60277188e+09 -9.37776409e+07 3.95191645e+09 -4.70829271e+09 -2.46177165e+09 -4.56039225e+09 -7.93992470e+08 -8.16655636e+07 3.76973782e+09 -9.08427126e+09 -6.09854749e+09 -4.65843749e+09 -1.48576388e+09 -4.68138797e+09 -1.72650565e+09 -4.69891777e+09 -1.97562736e+09 -4.71059934e+09 -2.23499480e+09 -4.71594846e+09 -2.50641944e+09 -4.42101039e+09 -4.54571675e+09 -1.14122677e+08 2.86889713e+09 -2.80084535e+09 -7.82099856e+08 -2.66504274e+09 4.00901070e+08 -2.61028150e+09 1.43488010e+08 -2.55249512e+09 -1.18210338e+08 -2.49081706e+09 -3.86734024e+08 -2.42437763e+09 -6.64519762e+08 -2.35229997e+09 -9.53923972e+08 -2.27369646e+09 -1.25724240e+09 -2.18766569e+09 -1.57672685e+09 -1.87549609e+09 -2.64203849e+09 -1.24029156e+08 1.67849325e+09 -1.19853017e+09 1.80523257e+09 -1.16531683e+09 2.57100807e+09 -1.05267517e+09 2.41485897e+09 -9.38374074e+08 2.25331883e+09 -8.21181766e+08 2.08454921e+09 -6.99912313e+08 1.90674438e+09 -5.73414618e+08 1.71811882e+09 -4.40563062e+08 1.51689660e+09 -3.00249526e+08 1.30130212e+09 -7.82038540e+07 9.75944932e+08 -1.18945133e+08 7.25503053e+08 -4.04371633e+08 4.67927759e+09 -5.32351458e+08 5.18169473e+09 -4.16936466e+08 5.20607423e+09 -3.00663837e+08 5.22470483e+09 -1.82352761e+08 5.23713199e+09 -6.08799664e+07 5.24282753e+09 6.48323612e+07 5.24119466e+09 1.95826797e+08 5.23157204e+09 3.33119359e+08 5.21323717e+09 3.77535807e+08 5.23078206e+09 -9.68745092e+07 8.47733475e+07 -4.88785356e+08 6.97978522e+09 -7.77795340e+08 7.27355108e+09 -7.13425775e+08 7.44961876e+09 -6.48339394e+08 7.61967992e+09 -5.81812139e+08 7.78445872e+09 -5.13150941e+08 7.94451622e+09 -1.04487613e+09 1.05963629e+10 -3.23580133e+08 5.84167160e+09 -3.11342575e+08 8.47624298e+09 -1.27896547e+08 8.56633154e+09 -3.82806143e+07 -1.78824005e+08 -1.54180399e+09 8.90988194e+09 -3.10386751e+09 1.85239713e+10 -5.18278028e+07 -1.55699211e+08 -1.41785043e+09 7.95733347e+09 -1.46879433e+09 8.13116314e+09 -1.49082351e+09 8.33360258e+09 -2.99285069e+09 1.74314455e+10 -1.99479614e+09 7.36938870e+09 -1.56407705e+09 7.46016093e+09 3.68354549e+07 4.07415404e+08 -1.31027067e+09 5.70558358e+09 -1.56265583e+09 5.58723398e+09 -1.69782019e+09 5.50601527e+09 -1.83214291e+09 5.41610728e+09 1.06828394e+07 -4.58989944e+08 -4.17026729e+09 1.09336758e+10 -2.24077780e+09 5.07972594e+09 -2.38207399e+09 4.94026067e+09 -2.52736520e+09 4.78442451e+09 -2.92340468e+09 5.10218614e+09 6.87199970e+07 1.43884752e+09 -1.44006182e+08 3.92849650e+09 -1.68535115e+08 3.57535399e+09 -2.72875883e+06 3.40289694e+09 -2.86357456e+08 3.34593339e+09 -4.05641488e+08 3.10554260e+09 -5.27471358e+08 2.85167461e+09 -6.52878367e+08 2.58183604e+09 -7.82851774e+08 2.29353442e+09 -9.18347962e+08 1.98426718e+09 -1.17615407e+09 1.76024301e+09 1.09094738e+08 3.14168906e+09 2.22323525e+09 1.25550782e+09 2.39612066e+09 2.52031041e+09 2.35028802e+09 2.26732426e+09 2.30227553e+09 2.00756350e+09 2.25129980e+09 1.73830353e+09 2.19658324e+09 1.45688947e+09 2.13734901e+09 1.16071594e+09 2.07281677e+09 8.47209774e+08 2.00219896e+09 5.13815104e+08 1.93748195e+09 -6.21796486e+08 1.22185104e+08 4.23982378e+09 4.79170234e+09 -2.90645628e+08 5.13349934e+09 1.80201523e+09 5.18855063e+09 1.62017225e+09 5.23928853e+09 1.43415355e+09 5.28554664e+09 1.24194318e+09 5.32708400e+09 1.04158943e+09 5.36359129e+09 8.31188043e+08 5.39469594e+09 6.08867943e+08 5.41996650e+09 3.72779011e+08 5.53682953e+09 -1.58775354e+09 1.19280159e+08 4.29218626e+09 6.87927676e+09 -1.13764953e+09 7.28697149e+09 1.26359915e+09 7.43499409e+09 1.19458116e+09 7.57696192e+09 1.12399279e+09 7.71333306e+09 1.05092621e+09 7.84441831e+09 9.74506271e+08 7.97039797e+09 8.93881901e+08 8.09133593e+09 8.08218979e+08 8.20719156e+09 7.16694229e+08 8.43570325e+09 -1.44294243e+09 1.03966674e+08 3.15325686e+09 7.71416119e+09 -1.44236268e+09 8.07064754e+09 5.36910690e+08 8.26834866e+09 5.58672990e+08 8.45927662e+09 5.80715473e+08 8.64428103e+09 6.03100179e+08 8.82402632e+09 6.25889693e+08 8.99901461e+09 6.49146676e+08 9.16960486e+09 6.72933560e+08 9.33602895e+09 6.97312376e+08 9.57702287e+09 -8.37498486e+08 1.49331315e+08 1.25953721e+09 6.84313865e+09 -1.25963073e+09 7.15056198e+09 -4.90100308e+08 7.33376081e+09 -4.08871368e+08 7.51102751e+09 -3.27064060e+08 7.68321646e+09 -2.43871278e+08 7.85101049e+09 -1.58524729e+08 8.01494267e+09 -7.02867477e+07 8.17541462e+09 2.15566224e+07 8.33271186e+09 1.17701511e+08 8.45040281e+09 -5.64007879e+08] using numpy matrix 0.9999420421616735 ```
mikerife commented 2 years ago

Hi @raphael-31 I'm very new to Numpy/Scipy and so while this is a shot in the dark, I've done some testing that seems to confirm. Been lost in the Numpy/Scipy documentation and have not found an definitive statement about this. Anyway the MAPDL stiffness matrix is (I'm guessing here but it seems to be by the description) symmetric. One of the symmetric modal solution types was used to solve correct? The .toarray and .todense create unsymmetric matrices but with zeros (or blanks which are treated as zeros) on the lower triangular part. If you create full matrices as shown in this example (K + K transpose - K diag) then the results match (in my test model):

https://mapdldocs.pyansys.com/examples/gallery_examples/01-apdlmath-examples/mapdl_vs_scipy.html#sphx-glr-examples-gallery-examples-01-apdlmath-examples-mapdl-vs-scipy-py

About 2/3rd of the way down the page the creation of the full matrices are shown.

My gut tells me that Numpy/Scipy have some way of flagging a matrix as symmetric...or maybe some modifier on the matrix multiplication etc functions. But I've not found anything yet. Mike

raphael-31 commented 2 years ago

Hi @raphael-31 I'm very new to Numpy/Scipy and so while this is a shot in the dark, I've done some testing that seems to confirm. Been lost in the Numpy/Scipy documentation and have not found an definitive statement about this. Anyway the MAPDL stiffness matrix is (I'm guessing here but it seems to be by the description) symmetric. One of the symmetric modal solution types was used to solve correct? The .toarray and .todense create unsymmetric matrices but with zeros (or blanks which are treated as zeros) on the lower triangular part. If you create full matrices as shown in this example (K + K transpose - K diag) then the results match (in my test model):

https://mapdldocs.pyansys.com/examples/gallery_examples/01-apdlmath-examples/mapdl_vs_scipy.html#sphx-glr-examples-gallery-examples-01-apdlmath-examples-mapdl-vs-scipy-py

About 2/3rd of the way down the page the creation of the full matrices are shown.

My gut tells me that Numpy/Scipy have some way of flagging a matrix as symmetric...or maybe some modifier on the matrix multiplication etc functions. But I've not found anything yet. Mike

Thank you Mike for your answer, but as I said, after I would analyse the result from an acoustic simulation, and for that type of result I don't think the stiff and the mass matrix are still symmetric.

mikerife commented 2 years ago

Hi @raphael-31 Please then provide either the model or describe the steps necessary to reproduce the same, or similar, model. Based on the file name the model does not appear to actually be an acoustic model: cylinder_no_water.dat implies the structural only part of the model. Mike

raphael-31 commented 2 years ago

Yes sure I attached it here. The first test I've made were indeed on a structural model (as it dosen't work as I want I don't try on the acoustic model). I would like to have the matrix (Mass, stiffness and damping) which verifies the eigenproblem in python with the firsts eigenvalues and eigenvectors to analyse it.

Thank you for your help, Raphaël.

Here is the acoustic model as .mac file :

Finish
/clear !clear the database

/PREP7 !preprocessor
etcontrol,set !automatically set KEYOPT settings

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
db_name = 'database'

!geometric variables 
average_radius = 1.0
width_structure = 0.1
outer_water_radius = 2.0

!divisions
radial_division = 10 
angular_division = 10

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!define element types
ET,1,FLUID29
!ET,2,FLUID129
!R, 1, outer_water_radius !real constant for the radius of the element 2
!ET, RAD, 2, outer_water_radius

ET,2,PLANE182

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!define material properties
!water
MP,DENS,1,1000, ! kg m^-3
MP,SONC,1,1430, ! m s^-1

!glass
!steel
MP,EX,2,70e9,   ! Pa
MP,PRXY,2,0.22, 
MP,DENS,2,2500, ! kg m^-3

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!geometry generation
!inner radius
CYL4,0,0, 0, 0, average_radius - width_structure, 90
CYL4,0,0, 0, 90, average_radius - width_structure, 180
CYL4,0,0, 0, 180, average_radius - width_structure, 270
CYL4,0,0, 0, 270, average_radius - width_structure, 360

!outer radius
CYL4,0,0, average_radius - width_structure, 0, average_radius + width_structure, 90
CYL4,0,0, average_radius - width_structure, 90, average_radius + width_structure, 180
CYL4,0,0, average_radius - width_structure, 180, average_radius + width_structure, 270
CYL4,0,0, average_radius - width_structure, 270, average_radius + width_structure, 360

!water outer radius
CYL4,0,0, average_radius + width_structure, 0, outer_water_radius, 90
CYL4,0,0, average_radius + width_structure, 90, outer_water_radius, 180
CYL4,0,0, average_radius + width_structure, 180, outer_water_radius, 270
CYL4,0,0, average_radius + width_structure, 270, outer_water_radius, 360

APTN,all !partition all volumes

NUMMRG,ALL, , , ,LOW !merge all volumes (merge common points between coincident volumes)
NUMCMP,ALL   !renumber volumes from 1 to x

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!MESHING controls
CSYS,1 ! activate cylindrical csys with Z as axis of rotation

!radial inner water division
LSEL,S,LOC,X, 0.0 + tol ,average_radius - width_structure - tol ! radial line selections
LESIZE,all, , ,radial_division, , , , ,1 ! set line divitions

!radial structure division
LSEL,S,LOC,X, average_radius - width_structure + tol ,average_radius + width_structure - tol ! radial line selections
LESIZE,all, , ,radial_division, , , , ,1 ! set line divitions

!radial outer water division
LSEL,S,LOC,X, average_radius + width_structure + tol ,outer_water_radius - tol ! radial line selections
LESIZE,all, , ,radial_division, , , , ,1 ! set line divitions

!angular division 
LSEL,S,LOC,X, average_radius - width_structure - tol ,average_radius - width_structure + tol
LSEL,A,LOC,X, average_radius + width_structure - tol ,average_radius + width_structure + tol 
LSEL,A,LOC,X, outer_water_radius - tol ,outer_water_radius + tol 
LESIZE,all, , ,angular_division, , , , ,1
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!Meshing
!steel structure meshing

! angular divisions
LSEL,S,LOC,X,inner_radius-tol,inner_radius+tol
LSEL,A,LOC,X,outer_radius-tol,outer_radius+tol
TYPE,2 !select element type  
MAT,2 !select material type
ASEL, S, LOC, X, average_radius - width_structure + tol, average_radius + width_structure - tol

AMESH, ALL !sweep steel bodies

!water meshing
TYPE,1 !select element type  
MAT,1 !select material type

ASEL,S,LOC,X,0 + tol, average_radius - width_structure - tol ! select outer parts (no shaft)
!ASEL,A,LOC,X, average_radius + width_structure + tol, outer_water_radius - tol

AMESH, ALL !sweep water bodies

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! create named components
ESEL,S,MAT,,1 ! select fluid elements (type 1 for fluid as defined above, type 2 for solid)
CM,fluid,elem ! create a component with all fluid elements
NSLE,S
CM,fluid_nodes,NODE ! create new component with all fluid nodes
ESEL,S,MAT,,2 ! select solid elements
CM,solid,elem ! create a component with all solid elements
NSLE,S
CM,solid_nodes,NODE ! create new component with all fluid nodes

ALLSEL,all

ALLSEL,all
!TYPE,2 !select element type  

!ESURF, ALL !sweep water bodies

! Boundary conditions
! On fixe un des noeuds
! csys, 0
! NSEL, S, LOC, X, average_radius + width_structure - tol, average_radius + width_structure + tol
! d,all,all !constrain all DOF at the selected nodes

!NSEL, S, LOC, X, average_radius - width_structure - tol, average_radius + width_structure + tol 
!d, all, UZ
nsel,all ! reselect all nodes
ddele,all,pres ! release the pressure dof for acoustics analysis
ecpchg !automatically convert uncoupled/coupled acoustic element for acoustic analysis

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Solution
/solu
antype,2 ! modal analysis
modopt,unsym,12,1.,1200.,AUTO !unsymetric solver, freq range 1Hz to 1200Hz
outres,erase
outres,all,none
outres,nsol,all
!mxpand,,,,yes,,no ! expand requested element results, but not write them to file.mode
dmpopt,esav,no
!emat est la matrice
dmpopt,emat,yes
dmpopt,full,yes
cdwrite,db,db_name,dat,,,,,BLOCKED !write database for the whole model (elements, nodes, etc.)
solve
FINISH
raphael-31 commented 2 years ago

Hi,

the model works ? And the extraction on python via numpy or scipy works ?

Have a good day, Raphaël.

germa89 commented 11 months ago

I gave it another look to this issue:

It seems the frequencies obtained are different:

MAPDL natural frequencies

array([101.3110303 , 101.32428665, 296.45490122, 296.70733178,
       556.97816386, 561.90897435, 571.5910606 , 574.59101973,
       732.93760406, 877.17729425])

Scipy natural frequencies

array([0.15964054, 0.15965804, 0.15966841, 0.15969207, 0.15969681,
       0.15973586, 0.15973672, 0.15974388, 0.15975051, 0.15976315])

There are big differences. Probably because of the symmetry of matrix?

Pinging @mikerife for another look

Code

from ansys.mapdl.core import launch_mapdl

mapdl = launch_mapdl()

cmd = """Finish
/clear !clear the database

/PREP7 !preprocessor
etcontrol,set !automatically set KEYOPT settings

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
db_name = 'database'

!geometric variables 
average_radius = 1.0
width_structure = 0.1
outer_water_radius = 2.0

!divisions
radial_division = 10 
angular_division = 10

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!define element types
ET,1,FLUID29
!ET,2,FLUID129
!R, 1, outer_water_radius !real constant for the radius of the element 2
!ET, RAD, 2, outer_water_radius

ET,2,PLANE182

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!define material properties
!water
MP,DENS,1,1000, ! kg m^-3
MP,SONC,1,1430, ! m s^-1

!glass
!steel
MP,EX,2,70e9,   ! Pa
MP,PRXY,2,0.22, 
MP,DENS,2,2500, ! kg m^-3

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!geometry generation
!inner radius
CYL4,0,0, 0, 0, average_radius - width_structure, 90
CYL4,0,0, 0, 90, average_radius - width_structure, 180
CYL4,0,0, 0, 180, average_radius - width_structure, 270
CYL4,0,0, 0, 270, average_radius - width_structure, 360

!outer radius
CYL4,0,0, average_radius - width_structure, 0, average_radius + width_structure, 90
CYL4,0,0, average_radius - width_structure, 90, average_radius + width_structure, 180
CYL4,0,0, average_radius - width_structure, 180, average_radius + width_structure, 270
CYL4,0,0, average_radius - width_structure, 270, average_radius + width_structure, 360

!water outer radius
CYL4,0,0, average_radius + width_structure, 0, outer_water_radius, 90
CYL4,0,0, average_radius + width_structure, 90, outer_water_radius, 180
CYL4,0,0, average_radius + width_structure, 180, outer_water_radius, 270
CYL4,0,0, average_radius + width_structure, 270, outer_water_radius, 360

APTN,all !partition all volumes

NUMMRG,ALL, , , ,LOW !merge all volumes (merge common points between coincident volumes)
NUMCMP,ALL   !renumber volumes from 1 to x

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!MESHING controls
CSYS,1 ! activate cylindrical csys with Z as axis of rotation

!radial inner water division
LSEL,S,LOC,X, 0.0 + tol ,average_radius - width_structure - tol ! radial line selections
LESIZE,all, , ,radial_division, , , , ,1 ! set line divitions

!radial structure division
LSEL,S,LOC,X, average_radius - width_structure + tol ,average_radius + width_structure - tol ! radial line selections
LESIZE,all, , ,radial_division, , , , ,1 ! set line divitions

!radial outer water division
LSEL,S,LOC,X, average_radius + width_structure + tol ,outer_water_radius - tol ! radial line selections
LESIZE,all, , ,radial_division, , , , ,1 ! set line divitions

!angular division 
LSEL,S,LOC,X, average_radius - width_structure - tol ,average_radius - width_structure + tol
LSEL,A,LOC,X, average_radius + width_structure - tol ,average_radius + width_structure + tol 
LSEL,A,LOC,X, outer_water_radius - tol ,outer_water_radius + tol 
LESIZE,all, , ,angular_division, , , , ,1
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!Meshing
!steel structure meshing

! angular divisions
LSEL,S,LOC,X,inner_radius-tol,inner_radius+tol
LSEL,A,LOC,X,outer_radius-tol,outer_radius+tol
TYPE,2 !select element type  
MAT,2 !select material type
ASEL, S, LOC, X, average_radius - width_structure + tol, average_radius + width_structure - tol

AMESH, ALL !sweep steel bodies

!water meshing
TYPE,1 !select element type  
MAT,1 !select material type

ASEL,S,LOC,X,0 + tol, average_radius - width_structure - tol ! select outer parts (no shaft)
!ASEL,A,LOC,X, average_radius + width_structure + tol, outer_water_radius - tol

AMESH, ALL !sweep water bodies

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! create named components
ESEL,S,MAT,,1 ! select fluid elements (type 1 for fluid as defined above, type 2 for solid)
CM,fluid,elem ! create a component with all fluid elements
NSLE,S
CM,fluid_nodes,NODE ! create new component with all fluid nodes
ESEL,S,MAT,,2 ! select solid elements
CM,solid,elem ! create a component with all solid elements
NSLE,S
CM,solid_nodes,NODE ! create new component with all fluid nodes

ALLSEL,all

ALLSEL,all
!TYPE,2 !select element type  

!ESURF, ALL !sweep water bodies

! Boundary conditions
! On fixe un des noeuds
! csys, 0
! NSEL, S, LOC, X, average_radius + width_structure - tol, average_radius + width_structure + tol
! d,all,all !constrain all DOF at the selected nodes

!NSEL, S, LOC, X, average_radius - width_structure - tol, average_radius + width_structure + tol 
!d, all, UZ
nsel,all ! reselect all nodes
ddele,all,pres ! release the pressure dof for acoustics analysis
ecpchg !automatically convert uncoupled/coupled acoustic element for acoustic analysis

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Solution
/solu
antype,2 ! modal analysis
modopt,unsym,12,1.,1200.,AUTO !unsymetric solver, freq range 1Hz to 1200Hz
outres,erase
outres,all,none
outres,nsol,all
!mxpand,,,,yes,,no ! expand requested element results, but not write them to file.mode
dmpopt,esav,no
!emat est la matrice
dmpopt,emat,yes
dmpopt,full,yes
cdwrite,db,db_name,dat,,,,,BLOCKED !write database for the whole model (elements, nodes, etc.)
solve
FINISH
"""
mapdl.clear()
out = mapdl.input_strings(cmd)

#
mm = mapdl.math

k = mm.stiff()
M = mm.mass()

nev = 10
A = mm.mat(k.nrow, nev)
ev = mm.eigs(nev, k, M, phi=A, fmin=1.0).asarray()

## Using scipy
pk = k.asarray()
pm = M.asarray()

# import numpy as np

# w, f = np.linalg.eig(np.linalg.inv(M)@k)
import scipy
import numpy as np
import math

pkd = scipy.sparse.diags(pk.diagonal())
pK = pk + pk.transpose() - pkd
pmd = scipy.sparse.diags(pm.diagonal())
pm = pm + pm.transpose() - pmd

from scipy.sparse.linalg import eigsh
vals, vecs = eigsh(A=pK, M=pm, k=nev, sigma=1, which="LA")

freqs = np.sqrt(vals) / (2 * math.pi)

assert np.allclose(ev, freqs), "The frequencies are different!"
mikerife commented 11 months ago

Hi @germa89 & @raphael-31 There was a FSI "surface force" flag missing from the MAPDL model. And SciPy needs the matrices cast to arrays (is cast the right word here?). I changed the model to be something we can explore in full i.e. 3 elements. SciPy does not return the eigenvalues in order so it's easier to show an example when you can show all the values.

from ansys.mapdl.core import launch_mapdl
import scipy
import numpy as np
import math
from scipy.sparse.linalg import eigs
import os
path = os.getcwd()

mapdl = launch_mapdl(run_location = path, additional_switches = '-smp', log_apdl= 'pymapdl_log.txt')

cmd = """
Finish
/clear 
/PREP7

ET,2,FLUID29
ET,3,FLUID29,,1,0  !press only
ET,1,PLANE182

!water
MP,DENS,2,1000, ! kg m^-3
MP,SONC,2,1430, ! m s^-1

!glass
!steel
MP,EX,1,70e9,    ! Pa
MP,PRXY,1,0.22, 
MP,DENS,1,2500,  ! kg m^-3

n,1,0,0,0
n,2,1,0,0
n,3,1,1,0
n,4,0,1,0

n,5,2,0,0
n,6,2,1,0

n,7,3,0,0
n,8,3,1,0

type,1
mat,1
real,1
secn,1
e,1,2,3,4

type,2
mat,2
real,2
secn,2
e,2,5,6,3

type,3
mat,2
real,3
secn,3
e,5,7,8,6

/solu

d,1,ux,0,,,,uy
d,4,ux,0,,,,uy

esel,s,elem,,2
nsel,s,node,,2,3
sf,all,fsi
allsel

antype,2 ! modal analysis
modopt,unsym,9,,,AUTO !unsymetric solver
outres,erase
outres,all,all
solve
FINISH
"""

mapdl.clear()
out = mapdl.input_strings(cmd)

This returns the following results: image

Then continuing on with SciPy:

mm = mapdl.math

k = mm.stiff()
M = mm.mass()

pk = k.asarray()
pm = M.asarray()

vals = eigs(A=pk.toarray(), M=pm.toarray(), k=9, sigma=1, which="LM", return_eigenvectors=False, tol = 0)
freqs = np.sqrt(vals) / (2 * math.pi)

freqs_abs = np.real(freqs)
print(np.sort(freqs_abs))
mapdl.exit()

The result of the print:

image

So no bug.
Mike

germa89 commented 11 months ago

I guess the missing flag is sf,all,fsi ??

germa89 commented 10 months ago

Closing because stale.