We have a user-subroutine that gets compiled and linked as a DLL to Finite Element Solver, Abaqus. Part of our subroutine is a matrix inversion subroutine that is inverting a 6x6 symmetric matrix. The subroutine we use is the one form https://people.sc.fsu.edu/~jburkardt/f_src/asa007/asa007.f90
The problem is that when the Abaqus simulation is run using a single cpu, the analysis runs fine. But when we use multiple cpus, it reports NANs.
Please note that the stiffness matrix being inverted can have a condition number of up to 1E7.
I wonder what can be the root cause of such problem. Can it be related to the matrix condition number?
>>The program is written in such a way that A and U can share storage.
Fortran specification: (Argument Association in Procedures)
Once an actual argument has been associated with a dummy argument, no action can be taken that affects the value or availability of the actual argument, except indirectly through the dummy argument.
Meaning if two (or more) dummy arguments represent overlapping memory, you are not permitted to change either.
Though in this case, since the single thread code works, I suspect this is not a case of vectorization of messing you up.
I suspect though that one of the callers uses a local "temp" array, and that that caller was compiled without -openmp (or recursive). Depending on version of Fortran, local arrays (non openmp, non recursive, non auto) have SAVE attribute-like behavior.