Two elements, both matrices, three modes

1. Structural model

Two prismatic frame elements meet at node 2. Nodes 1 and 3 are fully fixed, leaving only the three middle-node components out of nine global DOFs: u_x2, u_y2 and θ₂.

Exact inclined-element coordinates are (100/√2, 100/√2). The input file therefore uses 70.7106781186548, rather than the rounded 70.7: the latter would change the first element length by 0.0151%.

2. Executable model

The whole model fits into a short input file. The load is retained from the original example, but the method RunModal does not use it.

nodes
3
0,0,0
70.7106781186548,70.7106781186548,0
170.710678118655,70.7106781186548,0
elems_113
2
113,1,2,6,10000000,0.7,100
113,2,3,6,10000000,0.7,100
bcfix
2
1,1,0,0,0
1,3,0,0,0
bcforce_stat
1
10,2,0,-10,0

For a frame node, the DOF map is [ux, uy, thetaZ]. The free DOFs of node 2 therefore have global numbers [4 5 6], and the constrained DOFs are [1 2 3 7 8 9].

3. Local element matrices

Both elements have identical properties and length, hence identical local matrices. Numerical stiffness coefficients:

The local stiffness matrix, ordered as [u1 v1 θ1 u2 v2 θ2]:

k(local) = 10^6 ×
[ 0.600   0       0      -0.600   0       0
  0       0.012   0.600   0      -0.012   0.600
  0       0.600  40.000   0      -0.600  20.000
 -0.600   0       0       0.600   0       0
  0      -0.012  -0.600   0       0.012  -0.600
  0       0.600  20.000   0      -0.600  40.000 ]

For consistent mass, the factor preceding the dimensionless template is particularly convenient:

ρAL/420 = 0.7·6·100/420 = 1

m(local) =
[ 140      0       0      70      0       0
    0    156    2200       0     54   -1300
    0   2200   40000       0   1300  -30000
   70      0       0     140      0       0
    0     54    1300       0    156   -2200
    0  -1300  -30000       0  -2200   40000 ]

4. Rotation into global axes

For the first element, c=s=1/√2. The book factors out 1/√2; the implementation writes the same direction cosines directly:

T1 = 1/sqrt(2) ×
[ 1   1   0   0   0   0
 -1   1   0   0   0   0
  0   0  sqrt(2)  0   0   0
  0   0   0   1   1   0
  0   0   0  -1   1   0
  0   0   0   0   0  sqrt(2) ]

T2 = I

The same transformation applies to stiffness and mass:

k(global) = T^T k(local) T
m(global) = T^T m(local) T

This is where the inclined element couples global components x, y and rotation. The horizontal element needs no transformation.

5. Assembly in three free coordinates

Once the outer nodes are fixed, summing the free–free blocks of the two elements is sufficient. This is the same direct method that sums contributions by global DOF map in a larger model.

Element 1 contributions, 45° inclination

K1(ff) =
[ 306000   294000    424264.0687
  294000   306000   -424264.0687
  424264.0687 -424264.0687 40000000 ]

M1(ff) =
[ 148    -8     1555.6349
   -8   148    -1555.6349
 1555.6349 -1555.6349 40000 ]

Element 2 contributions, horizontal

K2(ff) =
[ 600000       0         0
       0   12000    600000
       0  600000  40000000 ]

M2(ff) =
[ 140    0      0
    0  156   2200
    0 2200  40000 ]

Reduced system matrices

Kff = K1(ff) + K2(ff)
    = 10^6 ×
[ 0.906000  0.294000  0.424264
  0.294000  0.318000  0.175736
  0.424264  0.175736 80.000000 ]

Mff = M1(ff) + M2(ff)
    =
[ 288        -8      1555.635
   -8       304       644.365
 1555.635   644.365 80000 ]

Rounding 0.424264 → 0.424, 0.175736 → 0.176, 1555.635 → 1556 and 644.365 → 644 gives exactly the matrices printed in the book.

6. Eigenvalues and frequencies

A nonzero solution exists when:

det(K_ff − ω²M_ff) = 0

Numerically, the generalized problem is solved: eig(Kff, Mff). Its eigenvalues equal λ=ω², after which we calculate ω=√λ and f=ω/(2π).

7. Mass-normalized modes

Each mode is scaled so that φ_i^TMφ_i=1. Columns of the modal matrix:

Phi =
[  0.02182975   0.00497889   0.05830737
  -0.05270167   0.00206232   0.02415170
   0            0.00340929  -0.00162915 ]

The printed book values [0.0218, -0.0527, 0], [0.00498, 0.00206, 0.00341] and [0.0583, 0.0241, -0.0016] are rounded versions of these columns.

Two required numerical checks:

Φ^TM_ffΦ = I

r_i = K_ffφ_i − λ_iM_ffφ_i ≈ 0

The eigenvector sign is arbitrary: φ and −φ describe the same mode. The test therefore first mass-normalizes each column, then aligns its sign with the reference before comparing components.

8. Reproducing in GNU Octave

addpath(pwd);
setup;

problem = StructFEProblem( ...
    fullfile('examples', 'cases', ...
        'CaseFig11.7p363 MarioPaz.txt'), ...
    struct('verbose', false, 'plotting', false));

model = problem.GetAnalysisModel();
result = problem.RunModal();

Kff = full(model.stiffness(result.freeDOFs, result.freeDOFs))
Mff = full(model.mass(result.freeDOFs, result.freeDOFs))
lambda = result.eigenvalues
omega = result.angularFrequenciesRadPerSec
frequencyHz = result.frequenciesHz
Phi = result.modeShapes(result.freeDOFs, :)

The automated check is in tests/test_mario_paz_modal.m and runs with the complete suite:

addpath(pwd); setup; run_octave_tests;

9. Practical experiments