PARDISO - trying to solve multiple matrices.

[potaman]

I am writing a non linear finite element code, where there is a time-stepping scheme and the matrix size and structure are different in each time step. I want to use PARDISO to solve the matrices.

Here is the calling scheme that I use.

do i=1,ntimesteps

call construct_matrix()
call pardiso with phase =11
call pardiso with phase =22

do till convergence
call pardiso with phase =33 !! different rhs. 



after convergence

call pardiso with phase = -1 and 0 !! free the memory L and U etc.

The problem is that during the second time step I get this error.

** Error in PARDISO  ( sequence_ido,parameters) error_num= 8
*** Input check: ia(neqns+1)_new 17425 _old 0 are incompatible
*** Input parameters: inconsistent error= 8 max_fac_store_in: 1
          matrix_number_in : 1 matrix_type_in   : 11
          ido_in           : 11 neqns_in         : 18288
          ia(neqns_in+1)-1 : 17425 nb_in            : 1

Any idea why this is a problem. The matrix structure and size changes at every time step , hence it is not possible to keep the same matrix alive. 

[Gennady_F_Intel]

the message "*** Input check: ia(neqns+1)_new 17425 _old 0 are incompatible" indicate thatits seem that you have some problem with ia array. Please check it accordingly with csr format.

[mecej4]

Not having seen the actual code, the best that one can do is to ask that you examine your source code using the information that PARDISO gave you: it remembers having been called with ia(neqns + 1) equal to zero. How did this come about?

[potaman]

The thing is that the matrix assembly is correct. i checked it by reading the matrices in to umfpack using scipy.

The main problem seems to be that the pardiso handle is not getting deleted. When I call the subroutine the first time , it seems to work okay. However on the second call, it seems to not have deallocated/deleted the previous handle. I can't post the entire code due to IP problems where I work.

[potaman]

I was trying to use the dss interface. While doing that, the problem seemed to be that a call to dss_delete did not seem to be doing anything...

[potaman]

Thanks guys for all the help. I figured out what was wrong. I did not know fortran 90 implicitly assumed variables initialized at declaration to be static. When I was calling my matrix construction the second time, it was not compressing the indices.

[mecej4]

> assumed variables initialized at declaration to be static

That you said "to be static" rather than "have the SAVE attribute" suggests a strong background in C, since initializations in C are executed at every entry to a function (or even a block) containing the declaration of the initialized variable.

Note, also, that in Fortran, despite the SAVE attribute, the variable need not have the value of the initialization expression when the containing subprogram is entered for the second, third, ... time. The value of the variable is whatever it was when control left the subprogram during the previous call.

[potaman]

I actually don't have a very strong background in any compiled language. I do most of my stuff in python. Thanks for all the help.