# NSASM Sparse Navier-Stokes Jacobian Assembly

NSASM is a FORTRAN90 library which carries out the finite element assembly of the jacobian matrix used in a Newton iteration to solve the steady state incompressible Navier Stokes equations in a 2D region, by Per-Olof Persson.

The finite element approximation uses the P2-P1 triangular element, also known as the Taylor-Hood element. The velocities are approximated quadratically, and the pressures linearly.

The sparse matrix format used is known as the compressed column format, in which the nonzero entries and their row indices are stored in order by columns.

The software assumes that a suitable mesh has been set up for the region, with arrays P and T defining that mesh, that an array defining the boundary constraints has been set up, and that an initial estimate of the flow has been made available.

The software was originally designed to be called, via MATLAB's mex facility, with the linear system solved by SUPER_LU. Thus, the Newton iteration could have the following form:

```        for ii = 1 : 8
[ K, L ] = nsasm ( p, t, np0, e, u, mu );
du = - lusolve ( K, L );
u = u + du;
disp ( sprintf ( '%20.10g  %20.10g', norm ( L, inf ), norm ( du, inf ) ) );
end
```

### Indexing Conventions

The software makes some particular assumptions about the numbering of nodes, variables, and local element nodes

(Indexing of nodes) It is assumed that the mesh was generated by starting with a set of nodes, triangulating them, and then computing the midsides of the triangles and adding these midsides to the list of nodes. In particular, the program therefore assumes that the numbering of the nodes reflects this process, so that all vertex nodes are numbered first, followed by the midside nodes.

(Indexing of global variables) It is assumed that every node has an associated horizontal velocity variable "U", that every node has an associated vertical velocity variable "V", and that only vertex nodes or "pressure nodes" have an associated pressure variable "P". The variables are stored in a single vector, first all the U's, then all the V's, then the P's. Thus variable 1 is the value of U at node 1, variable NP+1 is the value of V at node 1, and variable 2*NP+1 is the value of P at node 1.

Notice that there are even more global variables than you might think, since the boundary conditions and other constraints are handled by adding a Lagrange multiplier for each one. Thus the number of variables is actually 2*NP+NP0+NE.

(Numerical codes for U,V,P) In the E vector used to define constraints, the numeric codes 0, 1 and 2 are used for U, V and P variables respectively.

(Local numbering of nodes) When listing the 6 nodes that make up a triangle, the ordering should be as follows:

```        N3
|\
| \
|  \
N5   N4
|    \
|     \
|      \
N1--N6--N2
```

### Languages:

NSASM is available in a C version and a FORTRAN90 version and a MATLAB version.

### Related Data and Programs:

FEM2D_NAVIER_STOKES, a FORTRAN90 program which solves the 2D incompressible Navier Stokes equations in an arbitrary triangulated region. In order to run, it requires user-supplied routines that define problem data.

### Author:

Original C version by Per-Olof Persson;
FORTRAN90 version by John Burkardt.

### Reference:

Implementation of Finite Element-Based Navier-Stokes Solver,
April 2002.

### Examples and Tests:

The BIG test defines a relatively large sparse matrix, with NP = 2049 nodes, NT = 960 elements, NE = 1287 constraints, NP0 = 545, NU = 500.

The SMALL test defines a small sparse matrix, with NP = 25 nodes, NT = 8 elements, NE = 33 constraints, NP0 = 9, NU = 100.

### List of Routines:

• ASSEMBLE assembles the local stiffness and residual into global arrays.
• ASSEMBLE_CONSTR assembles the constraints.
• CH_CAP capitalizes a single character.
• CH_EQI is a case insensitive comparison of two characters for equality.
• CH_TO_DIGIT returns the integer value of a base 10 digit.
• FILE_COLUMN_COUNT counts the number of columns in the first line of a file.
• FILE_ROW_COUNT counts the number of row records in a file.
• GET_UNIT returns a free FORTRAN unit number.
• I4VEC_HEAP_D reorders an I4VEC into an descending heap.
• I4VEC_SORT_HEAP_A ascending sorts an I4VEC using heap sort.
• INIT_SHAPE evaluates the shape functions at the quadrature points.
• LOCALKL assembles the local stiffness matrix and residual.
• NSASM_INTERFACE is an interface for NSASM.