I'm computing the schur complement using iparm(36)=1. I do this in phases so that I first perform phase 331, then phase 332 in my own code, and then phase 333 in pardiso.
Prior to this, the schur complement matrix is calculated using phase 12.
It seems to me that one call (phase 331 then 333) with N right-hand sides is a lot faster than, say, 10 calls with N/10 right-hand sides. From the time usage it acutally seems like pardiso refactorizes the matrix for each subsequent call using phase 33x? Could this be?
I expected a single call (many rhs) to be faster than several calls (fewer rhs), but not by what seems to be a 10-fold (using the example above).
The documentation says that iparm(36)=1 should be used when only Schur complement needs to be computed. If you want to do solves afterwards, you should better use iparm(36)=2. So I'd say that your iparm parameters are a bit off for your purposes.
So, is there a reason why you want to do solves with iparm(36)=1? If there is one, let us know so that we dive a bit and check whether we truly don't save what we unexpectedly need to compute (factorization arrays) in case iparm(36)=1.
I hope this helps.