How to get Matrix^-1 * Matrix in Schur-complement

General Discussion
1

Hello members and developers of the FreeFEM：

While testing the Schur-complement method for solving a certain matrix(My thoughts are shown in the attachment), I noticed that FreeFEM++ does not support operations like Matrix^-1 * Matrix. so i modified the approach to use Matrix^-1 * real[int] instead, but this is undoubtedly very time-consuming—andthe numerical results confirm this. Is this method effectively impractical in FreeFEM?Additionally, I tried loading the “lapack” and “Schur-Complement” modules, but in both cases, the computation time was so longer than that of direct solution methods.

I would truly appreciate your guidance on this matter.

My program schematic:

Solve the following equation:
//[A B] [X1] = F ;
//[C D] [X2] = G ;
//========================
//Schur method
//========================
set(A, solver=sparsesolver);

//step(1) G-CA^-1F (G=0)
real[int] KinvF = Kfix^-1 * F;
real[int] CKinvF = C * KinvF;

//step(2) -CA^-1B+D (D=0)
matrix KinvB(B.n, B.m);for (int i = 0; i < B.m; i++) {
real[int] bi(B.n);
for (int j = 0; j < B.n; j++) {
bi[j] = B(j, i);
}
real[int] xi = Kfix^-1 * bi;

for (int j = 0; j < B.n; j++) {
    KinvB(j, i) = xi\[j\];  
}

}

matrix S = -1*C * KinvB;

//step(3) S*X = H
set(S, solver=sparsesolver);
real[int] X = S^-1 * CKinvF;

//===================
// time
//===================
Matrix A decompose time: 0.376S
step(1) calculate time: 0.107S
step(2) calculate time: 249.371S
step(3) calculate time: 107.208S

2

If you use PETSc, this can be done automatically for you, see, e.g., PCFieldSplitSetSchurFactType — PETSc 3.24.0 documentation.

3

Thank you for your response, Professor.
I have carefully read and studied the materials you provided, which has given me a better understanding of the subject.
However, since I have limited knowledge and practical experience with PETSc, I would like to further inquire: For my specific problem—using the A-PHI method to compute the electromagnetic field of an electric motor, while incorporating stator-rotor interface coupling, which requires solving at each time step—would this approach likely lead to a significant speedup compared to the direct method in FreeFEM (e.g., set(K, solver=sparsesolver);)?

4

You should definitely not compute the Schur complement yourself.