Reading and writing Nastran output4 files

This and other notebooks are available here: https://github.com/twmacro/pyyeti/tree/master/docs/tutorials.

First, do some imports:

import numpy as np
import os
from collections import OrderedDict
from pyyeti.nastran import op4

Generate a couple matrices to write to disk:

a = np.random.randn(4, 10)
b = np.random.randn(5, 14)

Writing op4 files

The nastran.op4 module can be used to read and write op4 files. We’ll demo writing first.

Using the dictionary option

If you don’t care about the order of the matrices in the file, you can use a regular dictionary as follows. By default, the file with be a native-endian binary and the matrices will be written in dense (non-sparse) format.

filename = 'rw_op4_demo.op4'
op4.write(filename, dict(a=a, b=b))

If you do care about order, you can specify the matrices using an OrderedDict. The following rewrites the file but ensures that “b” is before “a”:

dct = OrderedDict()
dct['b'] = b
dct['a'] = a
op4.write(filename, dct)

You can check the order with the op4.dir function:

op4.dir(filename);
b       ,      5 x 14    , form=2, mtype=2
a       ,      4 x 10    , form=2, mtype=2

If you need to specify the matrix form, you can include the desired value in the dictionary. Note, you only need to specify the forms you need. The default is None for each matrix, which means it will be automatically set. For example, if you want to mark matrix “a” as symmetric (which makes no sense here … “a” is not even square), but let “b” be automatically set:

dct = OrderedDict()
dct['b'] = b
dct['a'] = (a, 6)    # form=6 means "symmetric"
op4.write(filename, dct)
op4.dir(filename);
b       ,      5 x 14    , form=2, mtype=2
a       ,      4 x 10    , form=6, mtype=2

Using the list option

Providing the inputs via lists is an alternative method for writing op4 files. In fact, if you need to write multiple matrices of the same name, you must use the list option. For example, the following writes two matrices named “a” to the output file:

c = np.random.randn(3, 3)
c = c + c.T   # make it actually symmetric
names = ['a', 'b', 'a']
matrices = [a, b, c]
op4.write(filename, names, matrices)
op4.dir(filename);
a       ,      4 x 10    , form=2, mtype=2
b       ,      5 x 14    , form=2, mtype=2
a       ,      3 x 3     , form=6, mtype=2

The forms can also be input via a list. Unlike for dictionaries, if a form entry is needed for one matrix, an entry is needed for all matrices. Just use None’s where you want the automatic setting:

op4.write(filename, names, matrices, forms=[6, None, None])
op4.dir(filename);
a       ,      4 x 10    , form=6, mtype=2
b       ,      5 x 14    , form=2, mtype=2
a       ,      3 x 3     , form=6, mtype=2

Writing ASCII (text) op4 files:

As noted above, matrices are written in binary by default. To write in ASCII format, set binary to False:

dct = OrderedDict()
dct['a'] = a
dct['b'] = b
op4.write(filename, dct, binary=False)

Show the first 5 lines to see the ASCII format:

with open(filename) as f:
    for i in range(5):
        print(f.readline().strip())
10       4       2       2A       1P,3E23.16
1       1       4
1.1049687378816170E-01 6.1381448890413326E-01 5.1821525908682664E-01
1.4992467058688332E+00
2       1       4

When writing in ASCII, you can specify the number of digits (also marking “a” as symmetric for demonstration):

dct = OrderedDict()
dct['a'] = (a, 6)
dct['b'] = b
op4.write(filename, dct, binary=False, digits=9)
with open(filename) as f:
    for i in range(5):
        print(f.readline().strip())
10       4       6       2A       1P,5E16.9
1       1       4
1.104968738E-01 6.138144889E-01 5.182152591E-01 1.499246706E+00
2       1       4
-1.244264092E+00-4.746130818E-01-7.858781239E-01 2.148196359E-01

For comparison, rewrite the ASCII file using the “bigmat” sparse format:

op4.write(filename, dct, binary=False, digits=9, sparse='bigmat')
with open(filename) as f:
    for i in range(5):
        print(f.readline().strip())
10      -4       6       2A       1P,5E16.9
1       0      10
9       1
1.104968738E-01 6.138144889E-01 5.182152591E-01 1.499246706E+00
2       0      10

Rewrite again, but this time using the “nonbigmat” sparse format:

op4.write(filename, dct, binary=False, digits=9, sparse='nonbigmat')
with open(filename) as f:
    for i in range(5):
        print(f.readline().strip())
10       4       6       2A       1P,5E16.9
1       0       9
589825
1.104968738E-01 6.138144889E-01 5.182152591E-01 1.499246706E+00
2       0       9

Reading op4 files

Use op4.read to read op4 files.

Reading into an OrderedDict

op4.write(filename, dict(a=a, b=b))
dct = op4.read(filename)
assert np.all(dct['a'] == a)
assert np.all(dct['b'] == b)

The “pretty-print” class from pyYeti can help in viewing the dictionary:

from pyyeti.pp import PP
PP(dct);
<class 'collections.OrderedDict'>[n=2]
    'a': float64 ndarray 40 elems: (4, 10)
    'b': float64 ndarray 70 elems: (5, 14)

If you need to read in the form and type, set the justmatrix option to False (or use the load function … it is the same as read except the default on justmatrix). In this case, each dictionary entry is a 3-element tuple: (matrix, form, type):

dct = op4.read(filename, justmatrix=False)
assert np.all(dct['a'][0] == a)
assert np.all(dct['b'][0] == b)
PP(dct);
<class 'collections.OrderedDict'>[n=2]
    'a': [n=3]: (float64 ndarray: (4, 10), 2, 2)
    'b': [n=3]: (float64 ndarray: (5, 14), 2, 2)

Reading into lists

If the op4 file has duplicate names, you’ll have to read the variables into a list:

b2 = b + 10
op4.write(filename, ['b', 'a', 'b'], [b, a, b2])
names, mats, forms, types = op4.read(filename, into='list')
names
['b', 'a', 'b']
assert np.all(b == mats[0])
assert np.all(a == mats[1])
assert np.all(b2 == mats[2])

Reading and writing sparse matrices

There are two aspects to “sparse” for this module:

  1. There’s the on-disk format. Matrices can be written in “dense” or one of the two sparse formats: either “bigmat” or “nonbigmat”. The writing format is controlled by the sparse option of the write function.

  2. There’s also the in-memory format. Matrices can be read into regular numpy.ndarray matrices or into scipy.sparse sparse matrices. This is controlled by the sparse option of the read or load functions.

These two aspects are independent of each other. For example, numpy.ndarray matrices can be written in a sparse format and matrices written in the dense format can be read into scipy.sparse matrices.

To work with sparse matrices, we’ll need the scipy.sparse module:

import scipy.sparse as sp

First, create a 5 million by 5 million sparse matrix with just 3 elements to experiment with:

data = [2.3, 5, -100.4]
rows = [2, 500000, 350000]
cols = [3750000, 500000, 4999999]
a = sp.csr_matrix((data, (rows, cols)), shape=(5000000, 5000000))
a
<Compressed Sparse Row sparse matrix of dtype 'float64'
    with 3 stored elements and shape (5000000, 5000000)>
print(a)
<Compressed Sparse Row sparse matrix of dtype 'float64'
    with 3 stored elements and shape (5000000, 5000000)>
  Coords    Values
  (2, 3750000)      2.3
  (350000, 4999999) -100.4
  (500000, 500000)  5.0

Save matrix to op4 file. Note: sparse='bigmat' is the default for scipy.sparse matrices. But, it doesn’t hurt to specify it explicitly:

op4.write(filename, dict(a=a), sparse='bigmat')

Read matrix back in. The default when reading any matrix is to create regular numpy.ndarray matrices. Therefore, the sparse option is required here:

dct = op4.read(filename, sparse=True)
dct
OrderedDict([('a',
              <COOrdinate sparse matrix of dtype 'float64'
                    with 3 stored elements and shape (5000000, 5000000)>)])

All sparse matrices are returned in the COO format by default. To get the sparse matrix in the CSR format (for example) instead of the COO format, you can specify the sparse option as a two-tuple:

dct = op4.read(filename, sparse=(True, sp.coo_matrix.tocsr))
dct
OrderedDict([('a',
              <Compressed Sparse Row sparse matrix of dtype 'float64'
                    with 3 stored elements and shape (5000000, 5000000)>)])

Translate an op4 file to another op4 format

As a final example, these tools can be used to rewrite an op4 file in a different format very easily. This examples translates a sparse format binary op4 file to a simpler ascii format while preserving the matrix forms. For demonstration, we’ll define “a” as symmetric (form=6) even though it’s not even square:

a = np.random.randn(4, 10)
b = np.random.randn(5, 14)
dct = OrderedDict()
dct['b'] = b
dct['a'] = [a, 6]
op4.write(filename, dct, sparse='bigmat')

Translate it to simple non-sparse ascii, preserving the forms:

dct = op4.load(filename)
asciifile = filename.replace('.op4', '_ascii.op4')
op4.write(asciifile, dct, binary=False)

Check that the order and forms are the same:

op4.dir(filename);
b       ,      5 x 14    , form=2, mtype=2
a       ,      4 x 10    , form=6, mtype=2
op4.dir(asciifile);
b       ,      5 x 14    , form=2, mtype=2
a       ,      4 x 10    , form=6, mtype=2
os.remove(filename)
os.remove(asciifile)