Photo by James Owen on Unsplash A surprisingly popular blog-post written here is Exporting Stiffness Matrix from Ansys . A sensible follow up question is what can one do with the exported stiffness matrix? In a recent Xansys Forum post, a question was raised on how we can edit the stiffness matrix of a superelement and use it for our model. An approach presented below is to first create a superelement that has the same number of DOF and nodal location that will serve as a template. An APDL script can then be written to edit the stiffness matrix entries as desired before exporting to a new superelement *.SUB file for use in future models. The self-contained script below demonstrates this. /prep7 et ,1, 185 mp , ex, 1, 200e3 mp , prxy, 1, 0.33 w = 0.1 ! single element (note nodal locations) n , 1, w, -w, -w n , 2, w, w, -w n , 3, -w, w, -w n , 4, -w, -w, -w n , 5, w, -w, w n , 6, w, w, w n , 7, -w, w, w n , 8, -w, -w, w e , 1, 2, 3, 4, 5, 6, 7, 8 /solu antype , substr ! analy

For a simple 2D straight line, the equation is \( y = m \cdot x + c\). Writing it matrix form, we have:

\(y = \begin{bmatrix} x & 1 \end{bmatrix} \cdot \begin{bmatrix} m\\c \end{bmatrix}\)

With multiple data points in space (e.g. \((x_1,y_1), (x_2,y_2),\dots\) ) the equation looks like this...

\( \begin{bmatrix}y_1\\y_2 \\y_3\\ \vdots \end{bmatrix} = \begin{bmatrix} x_1 & 1\\ x_2 & 1 \\ x_3 & 1 \\ \vdots & \vdots \end{bmatrix} \cdot \begin{bmatrix} m\\c \end{bmatrix} \)

This matrix looks similar to \(y = A \cdot b\) where our

\(A = \begin{bmatrix} x_1 & 1 \\ x_2 & 1 \\ x_3 & 1 \\ \vdots & \vdots \end{bmatrix} \) and \(b = \begin{bmatrix} m\\c \end{bmatrix} \)

If \(A\) is a square matrix, we can inverse \(A\) to get \(A^{-1} \cdot y = b\). However when the \(A\) matrix is not square, it is over constrained. There are more equations than unknowns like the problem at hand above. That's where the Moore-Penrose inverse does it's magic. From Wikipedia, the pseudoinverse is written as:

\(A^{+} = (A^{T} \cdot A)^{-1} \cdot A^{T} \)

where \(A^{+} \cdot A = I \)

Rewriting the linear equation to solve for \(b\):

\((A^{T} \cdot A)^{-1} \cdot A^{T} \cdot y = b\)

To program it in Ansys efficiently, the equation is rearranged as:

\((A^{T} \cdot A) \cdot b\ = A^{T} \cdot y \)

**APDL script**

Here's an APDL script that works through it.

!! Known answers to shoot for

!! y = mm*x+cc

mm = 12

cc = 34

! Setup A Matrix and fake data with known answers

**/prep7**

nval = 20

*dim, Amat,, nval, 2

*dim, y,, nval, 1

*do,ct,1,nval

Amat(ct,1) = ct ! x variable

Amat(ct,2) = 1

y(ct,1) = (mm*ct+cc)*(1+0.01*rand(-1,1))

**! adds random noise**

*enddo

! Calculates Transpose

*dim, AmatT,, 2, nval

*mfun,AmatT(1,1),tran, Amat(1,1)

! Computes LHS Matrix

*dim, lhs,, 2,2

*moper, lhs(1,1), AmatT(1,1), mult, Amat(1,1) ! LHS = \((A^{T} \cdot A) \)

! Computes RHS Matrix

*dim, rhs,, nval,1

*moper, rhs(1,1), AmatT(1,1), mult, y(1,1)

**! RHS =**

**\(A^{T} \cdot y \)**

! Solves matrix

*dim, coef,, 2,1

*moper, coef(1,1), lhs(1,1), solv, rhs(1,1) ! LHS*coef = RHS, solve for coef

! Prints out Gradient & Intercept

*stat, coef

**Resources**

Here's the text file of the above script:

**[Link]**

Short lecture on it on YouTube:

**[Link]**

keep it up

ReplyDeleteThank you a lot, you solve my big problem

ReplyDelete